001: /*
002: * $RCSfile: NeuQuantOpImage.java,v $
003: *
004: * Copyright (c) 2005 Sun Microsystems, Inc. All rights reserved.
005: *
006: * Use is subject to license terms.
007: *
008: * $Revision: 1.2 $
009: * $Date: 2005/05/10 01:03:23 $
010: * $State: Exp $
011: */
012: package com.sun.media.jai.opimage;
013:
014: import java.awt.Image;
015: import java.awt.Rectangle;
016: import java.awt.image.DataBuffer;
017: import java.awt.image.Raster;
018: import java.awt.image.RenderedImage;
019: import java.util.ArrayList;
020: import java.util.Hashtable;
021: import java.util.LinkedList;
022: import java.util.ListIterator;
023: import java.util.Map;
024: import javax.media.jai.ImageLayout;
025: import javax.media.jai.LookupTableJAI;
026: import javax.media.jai.OpImage;
027: import javax.media.jai.PixelAccessor;
028: import javax.media.jai.PlanarImage;
029: import javax.media.jai.iterator.RandomIter;
030: import javax.media.jai.iterator.RandomIterFactory;
031: import javax.media.jai.ROI;
032: import javax.media.jai.ROIShape;
033: import javax.media.jai.UnpackedImageData;
034: import com.sun.media.jai.util.ImageUtil;
035:
036: /**
037: * An <code>OpImage</code> implementing the "ColorQuantizer" operation as
038: * described in <code>javax.media.jai.operator.ExtremaDescriptor</code>
039: * based on the median-cut algorithm.
040: *
041: * This is based on a java-version of Anthony Dekker's implementation of
042: * NeuQuant Neural-Net Quantization Algorithm
043: *
044: * NEUQUANT Neural-Net quantization algorithm by Anthony Dekker, 1994.
045: * See "Kohonen neural networks for optimal colour quantization"
046: * in "Network: Computation in Neural Systems" Vol. 5 (1994) pp 351-367.
047: * for a discussion of the algorithm.
048: *
049: * Any party obtaining a copy of these files from the author, directly or
050: * indirectly, is granted, free of charge, a full and unrestricted irrevocable,
051: * world-wide, paid up, royalty-free, nonexclusive right and license to deal
052: * in this software and documentation files (the "Software"), including without
053: * limitation the rights to use, copy, modify, merge, publish, distribute, sublicense,
054: * and/or sell copies of the Software, and to permit persons who receive
055: * copies from any such party to do so, with the only requirement being
056: * that this copyright notice remain intact.
057: *
058: * @see javax.media.jai.operator.ExtremaDescriptor
059: * @see ExtremaCRIF
060: */
061: public class NeuQuantOpImage extends ColorQuantizerOpImage {
062: /** four primes near 500 - assume no image has a length so large
063: * that it is divisible by all four primes
064: */
065: protected static final int prime1 = 499;
066: protected static final int prime2 = 491;
067: protected static final int prime3 = 487;
068: protected static final int prime4 = 503;
069:
070: /* minimum size for input image */
071: protected static final int minpicturebytes = (3 * prime4);
072:
073: /** The size of the histogram. */
074: private int ncycles;
075:
076: /* Program Skeleton
077: ----------------
078: [select samplefac in range 1..30]
079: [read image from input file]
080: pic = (unsigned char*) malloc(3*width*height);
081: initnet(pic,3*width*height,samplefac);
082: learn();
083: unbiasnet();
084: [write output image header, using writecolourmap(f)]
085: inxbuild();
086: write output image using inxsearch(b,g,r) */
087:
088: /* Network Definitions
089: ------------------- */
090:
091: private final int maxnetpos = maxColorNum - 1;
092: private final int netbiasshift = 4; /* bias for colour values */
093:
094: /* defs for freq and bias */
095: private final int intbiasshift = 16; /* bias for fractions */
096: private final int intbias = 1 << intbiasshift;
097: private final int gammashift = 10; /* gamma = 1024 */
098: private final int gamma = 1 << gammashift;
099: private final int betashift = 10;
100: private final int beta = intbias >> betashift; /* beta = 1/1024 */
101: private final int betagamma = intbias << (gammashift - betashift);
102:
103: /* defs for decreasing radius factor */
104: private final int initrad = maxColorNum >> 3;
105: private final int radiusbiasshift = 6; /* at 32.0 biased by 6 bits */
106: private final int radiusbias = 1 << radiusbiasshift;
107: private final int initradius = initrad * radiusbias; /* and decreases by a */
108: private final int radiusdec = 30; /* factor of 1/30 each cycle */
109:
110: /* defs for decreasing alpha factor */
111: private final int alphabiasshift = 10; /* alpha starts at 1.0 */
112: private final int initalpha = 1 << alphabiasshift;
113:
114: private int alphadec; /* biased by 10 bits */
115:
116: /* radbias and alpharadbias used for radpower calculation */
117: private final int radbiasshift = 8;
118: private final int radbias = 1 << radbiasshift;
119: private final int alpharadbshift = alphabiasshift + radbiasshift;
120: private final int alpharadbias = 1 << alpharadbshift;
121:
122: // typedef int pixel[4]; /* BGRc */
123: private int[][] network; /* the network itself - [maxColorNum][4] */
124:
125: private int[] netindex = new int[256]; /* for network lookup - really 256 */
126:
127: private int[] bias = new int[maxColorNum]; /* bias and freq arrays for learning */
128: private int[] freq = new int[maxColorNum];
129: private int[] radpower = new int[initrad]; /* radpower for precomputation */
130:
131: /**
132: * Constructs an <code>NeuQuantOpImage</code>.
133: *
134: * @param source The source image.
135: */
136: public NeuQuantOpImage(RenderedImage source, Map config,
137: ImageLayout layout, int maxColorNum, int upperBound,
138: ROI roi, int xPeriod, int yPeriod) {
139: super (source, config, layout, maxColorNum, roi, xPeriod,
140: yPeriod);
141:
142: colorMap = null;
143: this .ncycles = upperBound;
144: }
145:
146: protected synchronized void train() {
147:
148: // intialize the network
149: network = new int[maxColorNum][];
150: for (int i = 0; i < maxColorNum; i++) {
151: network[i] = new int[4];
152: int[] p = network[i];
153: p[0] = p[1] = p[2] = (i << (netbiasshift + 8))
154: / maxColorNum;
155: freq[i] = intbias / maxColorNum; /* 1/maxColorNum */
156: bias[i] = 0;
157: }
158:
159: PlanarImage source = getSourceImage(0);
160: Rectangle rect = source.getBounds();
161:
162: if (roi != null)
163: rect = roi.getBounds();
164:
165: RandomIter iterator = RandomIterFactory.create(source, rect);
166:
167: int samplefac = xPeriod * yPeriod;
168: int startX = rect.x / xPeriod;
169: int startY = rect.y / yPeriod;
170: int offsetX = rect.x % xPeriod;
171: int offsetY = rect.y % yPeriod;
172: int pixelsPerLine = (rect.width - 1) / xPeriod + 1;
173: int numSamples = pixelsPerLine
174: * ((rect.height - 1) / yPeriod + 1);
175:
176: if (numSamples < minpicturebytes)
177: samplefac = 1;
178:
179: alphadec = 30 + ((samplefac - 1) / 3);
180: int pix = 0;
181:
182: int delta = numSamples / ncycles;
183: int alpha = initalpha;
184: int radius = initradius;
185:
186: int rad = radius >> radiusbiasshift;
187: if (rad <= 1)
188: rad = 0;
189: for (int i = 0; i < rad; i++)
190: radpower[i] = alpha
191: * (((rad * rad - i * i) * radbias) / (rad * rad));
192:
193: int step;
194: if (numSamples < minpicturebytes)
195: step = 3;
196: else if ((numSamples % prime1) != 0)
197: step = 3 * prime1;
198: else {
199: if ((numSamples % prime2) != 0)
200: step = 3 * prime2;
201: else {
202: if ((numSamples % prime3) != 0)
203: step = 3 * prime3;
204: else
205: step = 3 * prime4;
206: }
207: }
208:
209: int[] pixel = new int[3];
210:
211: for (int i = 0; i < numSamples;) {
212: int y = (pix / pixelsPerLine + startY) * yPeriod + offsetY;
213: int x = (pix % pixelsPerLine + startX) * xPeriod + offsetX;
214:
215: try {
216: iterator.getPixel(x, y, pixel);
217: } catch (Exception e) {
218: continue;
219: }
220:
221: int b = pixel[2] << netbiasshift;
222: int g = pixel[1] << netbiasshift;
223: int r = pixel[0] << netbiasshift;
224:
225: int j = contest(b, g, r);
226:
227: altersingle(alpha, j, b, g, r);
228: if (rad != 0)
229: alterneigh(rad, j, b, g, r); /* alter neighbours */
230:
231: pix += step;
232: if (pix >= numSamples)
233: pix -= numSamples;
234:
235: i++;
236: if (i % delta == 0) {
237: alpha -= alpha / alphadec;
238: radius -= radius / radiusdec;
239: rad = radius >> radiusbiasshift;
240: if (rad <= 1)
241: rad = 0;
242: for (j = 0; j < rad; j++)
243: radpower[j] = alpha
244: * (((rad * rad - j * j) * radbias) / (rad * rad));
245: }
246: }
247:
248: unbiasnet();
249: inxbuild();
250: createLUT();
251: setProperty("LUT", colorMap);
252: setProperty("JAI.LookupTable", colorMap);
253: }
254:
255: private void createLUT() {
256: colorMap = new LookupTableJAI(new byte[3][maxColorNum]);
257: byte[][] map = colorMap.getByteData();
258: int[] index = new int[maxColorNum];
259: for (int i = 0; i < maxColorNum; i++)
260: index[network[i][3]] = i;
261: for (int i = 0; i < maxColorNum; i++) {
262: int j = index[i];
263: map[2][i] = (byte) (network[j][0]);
264: map[1][i] = (byte) (network[j][1]);
265: map[0][i] = (byte) (network[j][2]);
266: }
267: }
268:
269: /** Insertion sort of network and building of netindex[0..255]
270: * (to do after unbias)
271: */
272: private void inxbuild() {
273: int previouscol = 0;
274: int startpos = 0;
275: for (int i = 0; i < maxColorNum; i++) {
276: int[] p = network[i];
277: int smallpos = i;
278: int smallval = p[1]; /* index on g */
279: /* find smallest in i..maxColorNum-1 */
280: int j;
281: for (j = i + 1; j < maxColorNum; j++) {
282: int[] q = network[j];
283: if (q[1] < smallval) { /* index on g */
284: smallpos = j;
285: smallval = q[1]; /* index on g */
286: }
287: }
288: int[] q = network[smallpos];
289: /* swap p (i) and q (smallpos) entries */
290: if (i != smallpos) {
291: j = q[0];
292: q[0] = p[0];
293: p[0] = j;
294: j = q[1];
295: q[1] = p[1];
296: p[1] = j;
297: j = q[2];
298: q[2] = p[2];
299: p[2] = j;
300: j = q[3];
301: q[3] = p[3];
302: p[3] = j;
303: }
304: /* smallval entry is now in position i */
305: if (smallval != previouscol) {
306: netindex[previouscol] = (startpos + i) >> 1;
307: for (j = previouscol + 1; j < smallval; j++)
308: netindex[j] = i;
309: previouscol = smallval;
310: startpos = i;
311: }
312: }
313: netindex[previouscol] = (startpos + maxnetpos) >> 1;
314: for (int j = previouscol + 1; j < 256; j++)
315: netindex[j] = maxnetpos; /* really 256 */
316: }
317:
318: /** Search for BGR values 0..255 (after net is unbiased) and
319: * return colour index
320: */
321: protected byte findNearestEntry(int r, int g, int b) {
322: int bestd = 1000; /* biggest possible dist is 256*3 */
323: int best = -1;
324: int i = netindex[g]; /* index on g */
325: int j = i - 1; /* start at netindex[g] and work outwards */
326:
327: while (i < maxColorNum || j >= 0) {
328: if (i < maxColorNum) {
329: int[] p = network[i];
330: int dist = p[1] - g; /* inx key */
331: if (dist >= bestd)
332: i = maxColorNum; /* stop iter */
333: else {
334: i++;
335: if (dist < 0)
336: dist = -dist;
337: int a = p[0] - b;
338: if (a < 0)
339: a = -a;
340: dist += a;
341: if (dist < bestd) {
342: a = p[2] - r;
343: if (a < 0)
344: a = -a;
345: dist += a;
346: if (dist < bestd) {
347: bestd = dist;
348: best = p[3];
349: }
350: }
351: }
352: }
353:
354: if (j >= 0) {
355: int[] p = network[j];
356: int dist = g - p[1]; /* inx key - reverse dif */
357: if (dist >= bestd)
358: j = -1; /* stop iter */
359: else {
360: j--;
361: if (dist < 0)
362: dist = -dist;
363: int a = p[0] - b;
364: if (a < 0)
365: a = -a;
366: dist += a;
367: if (dist < bestd) {
368: a = p[2] - r;
369: if (a < 0)
370: a = -a;
371: dist += a;
372: if (dist < bestd) {
373: bestd = dist;
374: best = p[3];
375: }
376: }
377: }
378: }
379: }
380: return (byte) best;
381: }
382:
383: /** Unbias network to give byte values 0..255 and record
384: * position i to prepare for sort.
385: */
386: private void unbiasnet() {
387: for (int i = 0; i < maxColorNum; i++) {
388: network[i][0] >>= netbiasshift;
389: network[i][1] >>= netbiasshift;
390: network[i][2] >>= netbiasshift;
391: network[i][3] = i; /* record colour no */
392: }
393: }
394:
395: /** Move adjacent neurons by precomputed
396: * alpha*(1-((i-j)^2/[r]^2)) in radpower[|i-j|]
397: */
398:
399: private void alterneigh(int rad, int i, int b, int g, int r) {
400: int lo = i - rad;
401: if (lo < -1)
402: lo = -1;
403: int hi = i + rad;
404: if (hi > maxColorNum)
405: hi = maxColorNum;
406:
407: int j = i + 1;
408: int k = i - 1;
409: int m = 1;
410: while ((j < hi) || (k > lo)) {
411: int a = radpower[m++];
412: if (j < hi) {
413: int[] p = network[j++];
414: // try {
415: p[0] -= (a * (p[0] - b)) / alpharadbias;
416: p[1] -= (a * (p[1] - g)) / alpharadbias;
417: p[2] -= (a * (p[2] - r)) / alpharadbias;
418: // } catch (Exception e) {} // prevents 1.3 miscompilation
419: }
420: if (k > lo) {
421: int[] p = network[k--];
422: // try {
423: p[0] -= (a * (p[0] - b)) / alpharadbias;
424: p[1] -= (a * (p[1] - g)) / alpharadbias;
425: p[2] -= (a * (p[2] - r)) / alpharadbias;
426: // } catch (Exception e) {}
427: }
428: }
429: }
430:
431: /** Move neuron i towards biased (b,g,r) by factor alpha. */
432: private void altersingle(int alpha, int i, int b, int g, int r) {
433: /* alter hit neuron */
434: int[] n = network[i];
435: n[0] -= (alpha * (n[0] - b)) / initalpha;
436: n[1] -= (alpha * (n[1] - g)) / initalpha;
437: n[2] -= (alpha * (n[2] - r)) / initalpha;
438: }
439:
440: /** Search for biased BGR values. */
441: private int contest(int b, int g, int r) {
442: /* finds closest neuron (min dist) and updates freq */
443: /* finds best neuron (min dist-bias) and returns position */
444: /* for frequently chosen neurons, freq[i] is high and bias[i] is negative */
445: /* bias[i] = gamma*((1/maxColorNum)-freq[i]) */
446: int bestd = ~(((int) 1) << 31);
447: int bestbiasd = bestd;
448: int bestpos = -1;
449: int bestbiaspos = bestpos;
450:
451: for (int i = 0; i < maxColorNum; i++) {
452: int[] n = network[i];
453: int dist = n[0] - b;
454: if (dist < 0)
455: dist = -dist;
456: int a = n[1] - g;
457: if (a < 0)
458: a = -a;
459: dist += a;
460: a = n[2] - r;
461: if (a < 0)
462: a = -a;
463: dist += a;
464: if (dist < bestd) {
465: bestd = dist;
466: bestpos = i;
467: }
468: int biasdist = dist
469: - ((bias[i]) >> (intbiasshift - netbiasshift));
470: if (biasdist < bestbiasd) {
471: bestbiasd = biasdist;
472: bestbiaspos = i;
473: }
474: int betafreq = (freq[i] >> betashift);
475: freq[i] -= betafreq;
476: bias[i] += (betafreq << gammashift);
477: }
478: freq[bestpos] += beta;
479: bias[bestpos] -= betagamma;
480: return (bestbiaspos);
481: }
482: }
|