2.1.

Determining Seed Locations

Seed locations define the nuclei for growth of the spot function and the principal frequency of the screen. The seeds could be locations on a regular grid or irregular locations obtained from random perturbations of regular grids or from the application of an FM halftoning method to a uniform gray-level patch. Without loss of generality, we demonstrate the performance of the technique with seed locations obtained from a stochastic screen halftone of a uniform gray patch with area coverage $T$. Figure 2 shows a sample seed array.

Fig. 2

Example of seeds defined as activated pixels, obtained by applying a stochastic halftone threshold array to a patch with 6% area coverage.

2.2.

Determining Triangular Tessellation and Parameters

A triangular tessellation of a plane is a juxtaposition of triangular shapes in a pattern with no overlaps and no gaps. In order to tessellate the pixel area containing pixel locations, seeds are used as vertices of the triangles, and every seed is used as a vertex. If the seeds were arranged in a regular pattern, such as a hexagonal grid, then the tessellation would be regular and formed by congruent triangles. Figure 3 shows the triangular tessellation of a portion of the seed pattern from Fig. 2.

Fig. 3

Resulting triangular tessellation of a portion of the seed pattern from Fig. 2.

2.3.

Determining Triangle Membership

By definition, a triangular tessellation completely fills the space with no overlapping; therefore, every pixel in the image must belong to a triangle in the tessellation. For every pixel $P$ in the image, a determination is made as to which triangle encompasses it. Figure 4 shows an example of a pixel $P$ encompassed by triangle $P_1{P}_2{P}_3$, where $P_1$, $P_2$, and $P_3$ are the triangle vertices. There is no need to determine the encompassing triangle of a pixel if the pixel is a vertex. Once triangles are defined, relatively straightforward geometric relationships can be used to determine useful triangle parameters, such as triangle heights $H_1$, $H_2$, and $H_3$.

Fig. 4

Triangle $P_1{P}_2{P}_3$ encompassing pixel $P$. The triangle has heights $H_1$, $H_2$, and $H_3$ shown in dashed blue lines, and projected pixel heights $h_1$, $h_2$, and $h_3$ shown in dashed red lines.

2.4.

Determining Distances Between Pixel and Triangle Bases

Distances between pixels and triangle bases are determined from projections from the pixel location onto the triangle altitudes, or similarly, distances to the respective bases. An altitude of a triangle is a straight line through a vertex and perpendicular to the opposite side. The length of the altitude is the height with respect to the given vertex. The opposite side is called the base of the altitude. The distances $h_1$, $h_2$, and $h_3$ are the distances from $P$ to their respective bases, which are opposites $P_1$, $P_2$, and $P_3$, respectively, as illustrated in Fig. 4.

2.5.

Determining Spot Function Values

The distances calculated in Sec. 2.4 are used as inputs to a function that generates halftone thresholds. The function has parameters that can be adjusted to control its sharpness and slope. In our algorithm, the function is a weighted sum of three cosines that are functions of the three respective distances. Algebraic powers of the distances control sharpness of each dot touch point with its neighboring dots. Cosine weights control the dot touch sequencing, such that contact with neighboring dots can occur at different gray levels, thereby avoiding instability that may occur due to simultaneous touching. The spot function itself can be used to halftone an image, or a sampled version of it can be applied as a threshold array for efficient implementation in a printer. Given the pixel and triangle heights, we compute the spot function value $Q$ of that pixel according to

An additional set of parameters ${\gamma }_i>0$ can be introduced for control over the roundness of the dot sides and the sharpness of corner touch points.

Another level of control can be enabled by inverting the spot function. The inversion allows for well-controlled holes, which can be desirable for dark image subject matter.

Data normalization is an additional step that typically must be performed with spot functions. Once all pixels in the matrix are processed to obtain $Q$ values, the resulting thresholds are gray-level shifted and scaled to fit the data range, such as [0,255] for an 8-bit system or [0,1023] for a 10-bit image path. The effect of the growth control parameters is illustrated in Figs. 5Fig. 6Fig. 7Fig. 8Fig. 9 to 10, which show gray ramps halftoned with the aforementioned spot function having different parameter values. The images in the figures when printed at 2400 spots per inch have a principal frequency somewhat $>200$ cells per inch (cpi); from Fourier analysis of images, we determined it to be $\sim 230\mathrm{cpi}$. We note that although constructing the spot function is a computationally expensive process, it only has to be performed once for a given set of spot function parameters. The resulting fill order matrix can then be stored for future use and used within the framework of traditional screening methods to halftone arbitrary input images. We also note that, while halftone dot and hole symmetry is sometimes desired (e.g., a screen that transitions from dots with a predetermined shape to holes with the same shape at the midtones) in order to avoid undesirable textures or coarse patterns in the print,¹⁷ the stochastic nature of the proposed technique prevents it from resulting in perfectly symmetric patterns. For a technique that allows control of periodic and, if so desired, symmetrical halftone structures, we refer the reader to Sec. 3.

Fig. 5

Gray ramp halftoned with the spot function from Eq. (2) with $a_i=1$ and ${\gamma }_i=1$.

Fig. 6

Gray ramp halftoned with the spot function from Eq. (2) with $a_i=1$ and ${\gamma }_i=0.6$. Observe that ${\gamma }_i<1$ produces a pincushion appearance. This could improve stability in some marking processes.

Fig. 7

Gray ramp halftoned with the spot function from Eq. (2) with $a_i=1$ and ${\gamma }_i=1.5$. Observe that ${\gamma }_i>1$ tend to make dots rounder and delay touching at corners. This could be preferred for some marking processes or for subject matter consisting primarily of highlights.

Fig. 8

Gray ramp halftoned with the spot function from Eq. (2) when $a_i$ is varied with $H_i$, and ${\gamma }_i=1$. There are numerous ways in which the $a_i$ coefficients can be controlled. This image demonstrates increasing $a_i$ proportionally with $H_i$.

Fig. 9

Gray ramp halftoned with the spot function from Eq. (2) when $a_i$ is varied with $H_i$ and ${\gamma }_i=1$. This image demonstrates increasing $a_i$ inversely proportionally with $H_i$.

Fig. 10

Gray ramp halftoned with the spot function from Eq. (3) with $a_i=1$ and ${\gamma }_i=1$. Observe the well-controlled holes. This parameterization could be preferred for dark subject matter.

4.1.

Determining Dominant Orientation and Texture of Local Image Area

Image orientation and activity can be estimated in several ways. For example, the magnitude of the local gradient can serve as an estimate of the amount of spatial activity and its angle can be used to estimate the dominant orientation. Transform-based methods, such as principal component analysis, Gabor filtering, and Fourier, cosine, or wavelet decomposition, can also be used. In this paper, we propose the use of nonzero frequency (AC) coefficients in the discrete cosine transform (DCT) domain coefficients (which are readily available from JPEG-compressed data) to make that determination. Figure 16 shows a graphical representation of the first-level DCT basis functions used in JPEG, whose coefficients are used by the algorithm to determine dominant orientation and texture. Intuitively, a larger coefficient (0,1) indicates variation in the vertical direction of the local image block (typically $8\times 8\text{pixels}$). Similarly, the absolute values of coefficients (1,0) and (1,1) are proportional to the horizontal and diagonal edge activity of the local image block, respectively. Given these relationships, the magnitude of the AC coefficients can be used to estimate the principal orientation of activity within the image block. For simplicity, the current implementation of the algorithm approximates the local orientation to one of four possible angle values, from $-45$ to 90 deg in steps of 45 deg by using only the first-level JPEG decomposition coefficients, but a more accurate estimate of orientation could be achieved by using a larger number of coefficients. The estimated angle of orientation corresponds to that of the JPEG coefficient with the largest absolute value. If the largest coefficient is that in position (1,1), then the sign is also taken into account: a positive coefficient is indicative of a $-45\text{-}\mathrm{deg}$ angle, whereas a negative coefficient is indicative of a 45-deg orientation. More accurate estimates of local dominant directions can be achieved by using lower-level decomposition coefficients.

Fig. 16

Graphical representation of first-level DCT basis functions.

JPEG coefficients also carry information on the amount of local spatial variation. The current implementation of the algorithm estimates spatial activity by summing the first-level AC coefficients of the DCT decomposition of each image block (coefficients at coordinates (0,1), (1,0), and (1,1) in Fig. 16). Specifically, the first-level DCT coefficients are indicators of the main spatial component in each image block.

4.2.

Determining Seed Locations According to Orientation and Spatial Activity

A halftone screen is designed for each of a selected set of dominant image spatial activity orientations. Once the dominant direction of the spatial activity in the local image block is determined, its corresponding halftone screen is used to halftone an $8\times 8\text{pixel}$ block. The effects of this operation are twofold: first, the location of the seeds will be more aligned along the direction of the overall orientation of the image block; second, image blocks with larger spatial variation will have a larger number of seeds, effectively increasing the local principal frequency of the halftone in an orientation consistent with high image activity, thus improving detail rendition in those areas. For demonstration purposes, four anisotropic stochastic screens (one for each angle estimated in step 1) were designed, as illustrated in Fig. 17.

Fig. 17

Anisotropic stochastic screens used to determine seed locations (principal frequencies) for image blocks at (a) 0 deg, (b) 90 deg, (c) $-45\mathrm{deg}$, and (d) 45 deg.

The four anisotropic screens have a common subset of isotropically distributed dots to facilitate transitions between regions of different orientation; Fig. 18 shows this common subset of isotropic dots, which is also the identical output of the above four stochastic screens when halftoning a uniform input image at the low common gray level.

Fig. 18

Isotropic set of dots common to the four anisotropic screens depicted in Figs. 17(a) to 17(d).

Figure 19(b) shows the image resulting from halftoning the image of Fig. 19(a) with the algorithm described above.

Fig. 19

(a) Continuous tone image and (b) binary image resulting from halftoning the image in (a) using the algorithm of Sec. 4.

Figure 20 illustrates the enhanced edge and detail-rendering capabilities of the texture/orientation-adaptive seeding algorithm introduced in Sec. 4. Note how the number of seeds increases on edge regions, independent of the gray values of the image, effectively increasing the principal frequency of the halftone around areas of increased spatial activity. Also note the adaptive anisotropy, where the seed arrangement follows the edge orientation, even for cases where the edge has continuously changing orientation.

Fig. 20

Illustration of the enhanced edge and detail rendition capabilities of the algorithm from Sec. 4.

For comparison purposes, Fig. 21 shows the same halftone image area of Fig. 20(a) but obtained with the uniform isotropic principal frequency algorithm from Sec. 2. Figure 21(a) shows the halftoned image obtained with the common subset of seeds (Fig. 18), where the low-frequency rendition of the uniform sky would be printed in a stable low-noise manner, but the high-frequency content of the image is poorly represented by the low principal frequency of the screen. Figure 21(b) shows the opposite effect, where a union of the screens in Fig. 17, representing an isotropic high principal frequency screen, was used to better resolve fine image detail. The opposite problem occurs where the uniform sky is likely to print with more noise and less stability. These two screens represent the principal frequency range covered by the adaptive algorithm, which will produce the positive effects of these screens without incurring the negative [compare Fig. 20(a) to Figs. 21(a) and 21(b)].

Fig. 21

Performance of the traditional stochastically seeded clustered stochastic screen from Sec. 2 at (a) low-end and (b) high-end of the principal frequency range.

5.1.

Watermark Encoding

Encoding of a watermark is achieved by halftoning image content with two screens with different principal frequencies, the screens being selected by the watermarks. Figure 22 illustrates the steps involved in the creation of a watermarked halftone image, from a continuous-tone input image and a binary watermark mask. Recall that in the method introduced in Sec. 2, seeds were placed according to a binary patch corresponding to the stochastic screen halftoning of a uniform gray patch with area coverage $T$. In the watermark-embedding application, two different patches determine seed placement, the patches having been obtained from stochastic screen halftoning of two gray patches with area coverage $T_{\mathrm{W}}$ and $T_{\mathrm{NW}}$. At every pixel of the input image, seed location, and consequently, seed placement frequency (principal frequency), is determined based on the gray value of the pixel as well as on the presence or absence of watermark at that location. For locations where watermark (no watermark) is present, seeds are placed according to the binary patch obtained from the halftoning of a gray patch of value $T_{\mathrm{W}}(T_{\mathrm{NW}})$.

Fig. 22

Schematic illustration of the steps involved in the creation of a digital watermark.

Figures 23 and 24 provide an illustrative example of an encoded halftone image. Figure 23(a) shows the continuous-tone input image, and Fig. 23(b) shows a sample binary watermark mask. Figures 24(a) and 24(b) show pixel-size versions of the halftoned version of the image region highlighted in Fig. 23(a) for two different transition threshold gaps. In the case of Fig. 24(a), $T_{\mathrm{NW}}=6\%$ and $T_{\mathrm{W}}=15\%$. In contrast, Fig. 24(b) has transition thresholds $T_{\mathrm{NW}}=12\%$ and $T_{\mathrm{W}}=15\%$. It can be seen that, as the absolute difference between the transition thresholds $|T_{\mathrm{NW}}-T_{\mathrm{W}}|$ decreases, the perceptibility of the watermark decreases as well. The visibility of the watermark is thus easily controllable by adjusting the relative seed frequency between the watermark and nonwatermark areas.

Fig. 23

(a) Continuous tone image. (b) Binary watermarking mask.

Fig. 24

Binary images resulting from the halftoning of the image in Fig. 23(a) with the watermark in Fig. 23(b) using the introduced watermarking technique with (a) $T_{\mathrm{NW}}=6\%$, $T_{\mathrm{W}}=15\%$ and (b) $T_{\mathrm{NW}}=12\%$, $T_{\mathrm{W}}=15\%$.

5.2.

Watermark Decoding

Decoding of a watermark is achieved by identifying the principal frequency of the halftone screen used to halftone local regions in the halftoned image. In the context of this section, let us define the term dot as a set of 8-connected ON pixels and the term hole as a set of 8-connected OFF pixels. The decoding approach we propose is based on local counts of dots and holes. Consider the curves from Fig. 25(a), which illustrate the number of dots and holes for a traditional spot function designed according to the procedure from Sec. 2. For a spot function with transition threshold $T$, and as the gray level increases from 0, the number of dots increases steadily until gray level $T$ is reached, at which points it stabilizes. At a gray level where the dots start to touch, the number of dots decreases until all pixels in the screen are turned ON, at which point a single dot remains. In contrast, the number of holes remains steady at unity from a gray level of zero until a gray level where dots start to touch, at which point it increases steadily up to a level where dots grow so much that they start absorbing holes. From that point on, the number of holes decreases steadily until it reaches zero. Note that the shapes of the dot and hole count curves for a specific screen vary according to the choice of $T$ and the spot function utilized.

Fig. 25

Dot and hole count curves for (a) a single spot function with transition threshold $T$; and for two spot functions with thresholds $T_1$ and $T_2$ (b) when no erosion or dilation is applied and (c) after application of erosion and dilation operations.

Consider the case where two spot functions with two different transition threshold values $T_1$ and $T_2$ are used; the dot and hole count curves in this scenario are plotted in Fig. 25(b). Inspection of Fig. 25(b) reveals that the number of dots and holes at a portion of the midtones and shadows, respectively, can be used to discriminate between the screens, while the light shadows (levels below the smallest transition threshold, in this case $T_1$) and very dark shadows are not useable for decoding purposes. The useable gray-scale range can be increased by breaking the touch points between dots in the light shadows via morphological erosion and by filling in holes in the shadows via morphological dilation. The decision of whether to use erosion and dilation or not is made based on estimated area coverage information or on a comparison between local dot and hole count and predetermined thresholds. Figure 25(c) shows the new dot and hole count curves resulting from the application of the erosion and dilation operations.

Figure 26 shows the image resulting from the decoding of the invisible watermark shown in the image from Fig. 24(b). White pixels denote data identified as watermark. The retrieval of the watermark is not perfect because some portions of it fall in the unusable gray-scale range.

Fig. 26

Binary image resulting from the decoding of the embedded watermark from Fig. 24(b).

6.1.

Intuitive Description of the Successive Fill Strategy

The fill order for a scalar halftone spot function determined according to the procedure outlined in Sec. 2 establishes a threshold $t$ associated with each location within a triangle, which can be uniquely described by the distances from the specific location onto each of the triangle sides, namely $h_1$, $h_2$, and $h_3$. In traditional successive fill algorithms, the thresholds in the vector are determined by a pair of rectilinear coordinates within a rectangular cell based on simple modulo indexing into the cell; in contrast, the present vector elements are associated with the aforementioned orthogonal distances to sides of spatially varying triangles that are determined from two-dimensional point coordinates in space, as illustrated in Fig. 27.

Fig. 27

Vector of thresholds (in ascending order) associated with orthogonal distances between the location of the point and the sides of spatially varying triangles.

For a given input pixel, the successive fill algorithm is implemented as follows:

A filling strategy implemented in this manner results in a first colorant dot being grown from the seed locations, or triangle vertices, a second colorant dot being grown from the triangle center location, and a third coloration dot being grown from the midpoints of each of the triangle sides, as illustrated in Fig. 28.

Fig. 28

Origin of growth regions for image areas with color content in at least three different channels. The particular colorants in the figure are only used as an example.

6.2.

Algorithmic Implementation of the Successive Fill Strategy

The algorithmic implementation of the algorithm outlined above is described next. For a given input pixel having colorant values $\mathbf{f}=(f_1,f_2,f_3)$, we determine halftoned pixel values $\mathbf{b}=(b_1,b_2,b_3)$ by applying the following procedure:

Application of rule 3 to the lightest colorant is optional because the particular spot functions that we provided above will have the least compact form for this rule and, therefore, could print less reliably in this scenario. Figure 29 shows a schematic of this successive fill algorithm, which is significantly different from successive fill methods used for traditional stochastic screens.²⁷

Fig. 29

Schematic showing the threshold order used in the presented successive fill method.

Figures 30 and 31 show two- and three-colorant sweeps, respectively, halftoned with the algorithm described above. The dot-off-dot structures are maintained up to 50 and 33% area coverage, respectively. This dot-off-dot structure produces smoother textures and a larger gamut while using less ink compared to scalar halftoning techniques applied independently to each color channel. At the same time, and given the clustered nature of the screens, it offers improved stability. One potential problem associated with dot-off-dot screens is their sensitivity to color shift due to misregistration between the color planes in the print engine. Some degree of level blending can be built into the spot functions to alleviate this sensitivity. The specific implications of misregistration between color planes on the quality of the presented dot-off-dot halfotning technique are beyond the scope of the present paper. For an in-depth study on the effects of misregistration relative to induced color shifts in clustered color halftones, we refer the reader to Ref. 28.

Fig. 30

Two-colorant sweeps halftoned with the proposed algorithm. The dot-off-dot structures are maintained up to 50% area coverage per colorant.

Fig. 31

Three-colorant sweeps halftoned with the proposed algorithm. The dot-off-dot structures are maintained up to 33% area coverage per colorant.