001: package JSci.maths.wavelet;
002:
003: import JSci.maths.wavelet.*;
004: import JSci.maths.*;
005:
006: /****************************************************
007: * This class is used to encapsulate wavelet coefficients.
008: * @author Daniel Lemire
009: *****************************************/
010: public final class FWTCoef extends Object implements
011: NumericalConstants, Cloneable {
012: protected double[][] coefs;
013: final static double normalisation = 1.0 / SQRT2;
014:
015: public FWTCoef() {
016: }
017:
018: /**********************************************
019: ***********************************************/
020: public FWTCoef(double[][] v) {
021: coefs = v;
022:
023: }
024:
025: /********************************************
026: * Return a copy of this object
027: *********************************************/
028: public Object clone() {
029: try {
030: FWTCoef fwt = (FWTCoef) super .clone();
031: if (coefs != null)
032: fwt.coefs = ArrayMath.copy(coefs);
033: return (fwt);
034: } catch (CloneNotSupportedException cnse) {
035: throw new InternalError();
036: }
037: }
038:
039: /*************************************************
040: **************************************************/
041: public int getJ() {
042: return (coefs.length);
043: }
044:
045: /*************************************************
046: **************************************************/
047: public int dimension(int i) {
048: if ((i < 0) || (i >= coefs.length)) {
049: throw new IllegalArgumentException(
050: "This dimension doesn't exist : " + i + ", "
051: + coefs.length);
052: }
053: return (coefs[i].length);
054: }
055:
056: /*******************************************************
057: ********************************************************/
058: public double[][] getCoefs() {
059: return (coefs);
060: }
061:
062: /***************************
063: * Compute the L2 norm of the
064: * coefficients
065: ****************************/
066: public double[] norm() {
067: double[] ans = new double[coefs.length];
068: for (int j = 0; j < coefs.length; j++) {
069: ans[j] = ArrayMath.norm(coefs[j]);
070: }
071: return (ans);
072: }
073:
074: /***************************
075: * Compute the L2 norm of the
076: * coefficients at "scale" i.
077: * Wavelet coefficients are
078: * into the "scale" 1 to ... and
079: * the scale 0 is the coarsest
080: * scale containing scaling
081: * functions coefficients
082: ****************************/
083: public double norm(int i) {
084: if ((i < 0) || (i >= coefs.length)) {
085: throw new IllegalArgumentException("The integer parameter "
086: + i + " should be between 0 and "
087: + (coefs.length - 1));
088: }
089: double ans = ArrayMath.norm(coefs[i]);
090: return (ans);
091: }
092:
093: /************************************
094: * Compute the sum of the squares of
095: * the coefficients
096: *************************************/
097: private double[] sumSquares() {
098: double[] ans = new double[coefs.length];
099: for (int j = 0; j < coefs.length; j++) {
100: ans[j] = ArrayMath.sumSquares(coefs[j]);
101: }
102: return (ans);
103: }
104:
105: /************************************
106: * Compute the sum of the squares of
107: * the coefficients
108: *************************************/
109: public double sumSquares(int i) {
110: if ((i < 0) || (i >= coefs.length)) {
111: throw new IllegalArgumentException("The integer parameter "
112: + i + " should be between 0 and "
113: + (coefs.length - 1));
114: }
115: double ans = ArrayMath.sumSquares(coefs[i]);
116: return (ans);
117: }
118:
119: /************************************
120: *************************************/
121: public double mass(int i) {
122: if ((i < 0) || (i >= coefs.length)) {
123: throw new IllegalArgumentException("The integer parameter "
124: + i + " should be between 0 and "
125: + (coefs.length - 1));
126: }
127: double ans = ArrayMath.mass(coefs[i]);
128: return (ans);
129: }
130:
131: /************************************
132: *************************************/
133: private double[] variance() {
134: double[] ans = new double[coefs.length];
135: for (int j = 0; j < coefs.length; j++) {
136: ans[j] = ArrayMath.variance(coefs[j]);
137: }
138: return (ans);
139: }
140:
141: /************************************
142: *************************************/
143: public double variance(int i) {
144: if ((i < 0) || (i >= coefs.length)) {
145: throw new IllegalArgumentException("The integer parameter "
146: + i + " should be between 0 and "
147: + (coefs.length - 1));
148: }
149: double ans = ArrayMath.variance(coefs[i]);
150: return (ans);
151: }
152:
153: /**********************************************
154: ***********************************************/
155: public double sumEnergies() {
156: if (coefs.length <= 1) {
157: throw new IllegalArgumentException(
158: "No wavelet coefficients!");
159: }
160: double[] energies = sumSquares();
161: double ans = 0;
162: for (int k = 1; k < energies.length; k++) {
163: ans += energies[k];
164: }
165: return (ans);
166: }
167:
168: /******************************************************
169: *******************************************************/
170: public double entropy() {
171: if (coefs.length <= 1) {
172: throw new IllegalArgumentException(
173: "No wavelet coefficients!");
174: }
175: double se = sumEnergies();
176: int nombreDeCoefficients = 0;
177: for (int k = 1; k < coefs.length; k++) {
178: nombreDeCoefficients += coefs[k].length;
179: }
180: double[] er = new double[nombreDeCoefficients];
181: int pos = 0;
182: for (int k = 1; k < coefs.length; k++) {
183: for (int l = 0; l < coefs[k].length; l++) {
184: er[pos] = coefs[k][l] * coefs[k][l] / se;
185: pos++;
186: }
187: }
188: return (EngineerMath.icf(er));
189: }
190:
191: /**********************************************
192: ***********************************************/
193: public double sumVariance() {
194: if (coefs.length <= 1) {
195: throw new IllegalArgumentException(
196: "No wavelet coefficients!");
197: }
198: double[] variances = variance();
199: double ans = 0;
200: for (int k = 1; k < variances.length; k++) {
201: ans += variances[k];
202: }
203: return (ans);
204: }
205:
206: /***********************************************
207: ************************************************/
208: public double energyRatio(int i) {
209: if (coefs.length <= 1) {
210: throw new IllegalArgumentException(
211: "No wavelet coefficients!");
212: }
213: if ((i < 1) || (i >= coefs.length)) {
214: throw new IllegalArgumentException("The integer parameter "
215: + i + " should be between 0 and "
216: + (coefs.length - 1));
217: }
218: if (sumEnergies() == 0) {
219: if (coefs.length != 0) {
220: return (1 / coefs.length);
221: } else {
222: throw new IllegalArgumentException("No energy!");
223: }
224: }
225: return (sumSquares(i) / sumEnergies());
226: }
227:
228: /***********************************************
229: ************************************************/
230: public double varianceRatio(int i) {
231: if (coefs.length <= 1) {
232: throw new IllegalArgumentException(
233: "No wavelet coefficients!");
234: }
235: if ((i < 1) || (i >= coefs.length)) {
236: throw new IllegalArgumentException("The integer parameter "
237: + i + " should be between 0 and "
238: + (coefs.length - 1));
239: }
240: if (sumVariance() == 0) {
241: if (coefs.length != 0) {
242: return (1 / coefs.length);
243: } else {
244: throw new IllegalArgumentException("No energy!");
245: }
246: }
247: return (variance(i) / sumVariance());
248: }
249:
250: /***************************************************
251: * Compute the Shannon entropy.
252: ****************************************************/
253: public double icf() {
254: if (coefs.length <= 1) {
255: throw new IllegalArgumentException(
256: "No wavelet coefficients!");
257: }
258: double[] pe = new double[coefs.length - 1];
259: for (int j = 1; j < coefs.length; j++) {
260: pe[j - 1] = energyRatio(j);
261: }
262: return (EngineerMath.icf(pe));
263: }
264:
265: /***************************************************
266: ****************************************************/
267: public double varianceICF() {
268: if (coefs.length <= 1) {
269: throw new IllegalArgumentException(
270: "No wavelet coefficients!");
271: }
272: double[] pv = new double[coefs.length - 1];
273: for (int j = 1; j < coefs.length; j++) {
274: pv[j - 1] = varianceRatio(j);
275: }
276: return (EngineerMath.icf(pv));
277: }
278:
279: /***************************
280: ****************************/
281: public void setCoefs(double[][] v) {
282: coefs = v;
283: }
284:
285: /***************************
286: ****************************/
287: public void setCoefs(double[] v, int i) {
288: if ((i < 0) || (i >= coefs.length)) {
289: throw new IllegalArgumentException("The integer parameter "
290: + i + " should be between 0 and "
291: + (coefs.length - 1));
292: }
293: coefs[i] = v;
294: }
295:
296: /**********************************************
297: ***********************************************/
298: public void synthesize(Filter filtreprimaire, double[] param) {
299: if (coefs.length <= 1) {
300: throw new IllegalArgumentException(
301: "No synthesis possible : " + coefs.length);
302: }
303: double[] V0 = filtreprimaire.lowpass(coefs[0], param);
304: double[] W0 = filtreprimaire.highpass(coefs[coefs.length - 1],
305: param);
306: V0 = ArrayMath.scalarMultiply(normalisation, V0);
307: if (V0.length != W0.length) {
308: throw new IllegalArgumentException(
309: "Synthesis impossible : bad data/multiresolution?"
310: + coefs[0].length + ", "
311: + coefs[coefs.length - 1].length + ", "
312: + V0.length + ", " + W0.length);
313: }
314: V0 = ArrayMath.add(V0, W0);
315: double[][] c = new double[coefs.length - 1][];
316: for (int j = 1; j < coefs.length - 1; j++) {
317: c[j] = coefs[j];
318: }
319: c[0] = V0;
320: coefs = c;
321: }
322:
323: /**********************************************
324: ***********************************************/
325: public void synthesize(Filter filtreprimaire, double[] param,
326: int jmax) {
327: if ((jmax < 0) || (jmax > coefs.length - 1)) {
328: throw new IllegalArgumentException("The integer parameter "
329: + jmax + " must be between 0 and "
330: + (coefs.length - 1));
331: }
332: for (int j = 0; j < jmax; j++) {
333: synthesize(filtreprimaire, param);
334: }
335: }
336:
337: /**********************************************
338: ***********************************************/
339: public void synthesizeAll(Filter filtreprimaire, double[] param) {
340: synthesize(filtreprimaire, param, coefs.length - 1);
341: }
342:
343: /**********************************************
344: ***********************************************/
345: public void synthesize(Filter filtreprimaire) {
346: if (coefs.length <= 1) {
347: throw new IllegalArgumentException(
348: "No synthesis possible : " + coefs.length);
349: }
350: double[] V0 = filtreprimaire.lowpass(coefs[0]);
351: double[] W0 = filtreprimaire.highpass(coefs[coefs.length - 1]);
352: V0 = ArrayMath.scalarMultiply(normalisation, V0);
353: if (V0.length != W0.length) {
354: throw new IllegalArgumentException(
355: "Synthesis impossible : bad data/multiresolution?"
356: + coefs[0].length + ", "
357: + coefs[coefs.length - 1].length + ", "
358: + V0.length + ", " + W0.length);
359: }
360: V0 = ArrayMath.add(V0, W0);
361: double[][] c = new double[coefs.length - 1][];
362: for (int j = 1; j < coefs.length - 1; j++) {
363: c[j] = coefs[j];
364: }
365: c[0] = V0;
366: coefs = c;
367: }
368:
369: /**********************************************
370: ***********************************************/
371: public void synthesize(Filter filtreprimaire, int jmax) {
372: if ((jmax < 0) || (jmax > coefs.length - 1)) {
373: throw new IllegalArgumentException("The integer parameter "
374: + jmax + " must be between 0 and "
375: + (coefs.length - 1));
376: }
377: for (int j = 0; j < jmax; j++) {
378: synthesize(filtreprimaire);
379: }
380: }
381:
382: /**********************************************
383: ***********************************************/
384: public void synthesizeAll(Filter filtreprimaire) {
385: synthesize(filtreprimaire, coefs.length - 1);
386: }
387:
388: /**************************************************
389: ***************************************************/
390: public Signal rebuildSignal(Filter filtreprimaire) {
391: FWTCoef fwt = new FWTCoef(coefs);// copie
392: fwt.synthesizeAll(filtreprimaire);
393: return (new Signal(fwt.getCoefs()[0]));
394: }
395:
396: /**************************************************
397: ***************************************************/
398: public Signal rebuildSignal(Filter filtreprimaire, double[] param) {
399: FWTCoef fwt = new FWTCoef(coefs);// copie
400: fwt.synthesizeAll(filtreprimaire, param);
401: return (new Signal(fwt.getCoefs()[0]));
402: }
403:
404: /*********************************************
405: * Denoises by zero-ing any value above a given percentile cut-off.
406: * @param p percentile cut-off, must be between 0 and 1.
407: **********************************************/
408: public void denoise(double p) {
409: for (int k = 1; k < coefs.length; k++) {
410: coefs[k] = denoise(coefs[k], p);
411: }
412: }
413:
414: /*********************************************
415: * Denoises by zero-ing any value above a given percentile cut-off.
416: * @param p percentile cut-off, must be between 0 and 1.
417: * @param k the index of the coefficient array to denoise.
418: **********************************************/
419: public void denoise(double p, int k) {
420: coefs[k] = denoise(coefs[k], p);
421: }
422:
423: /**************************************
424: * Denoises by zero-ing any value above a given percentile cut-off.
425: * @param v an array to denoise.
426: * @param p percentile cut-off, must be between 0 and 1.
427: ***************************************/
428: public static double[] denoise(double[] v, double p) {
429: if (p == 0)
430: return (v);
431: double[] ans = v;
432: double seuil = ArrayMath.percentile(ArrayMath.abs(ans), 1 - p);
433: for (int k = 0; k < ans.length; k++) {
434: if (Math.abs(ans[k]) >= seuil) {
435: ans[k] = 0;
436: }
437: }
438: return (ans);
439: }
440:
441: /*********************************************
442: * Compresses by zero-ing any value below a given percentile cut-off.
443: * @param p percentile cut-off, must be between 0 and 1.
444: **********************************************/
445: public void compress(double p) {
446: for (int k = 1; k < coefs.length; k++) {
447: coefs[k] = compress(coefs[k], p);
448: }
449: }
450:
451: /*********************************************
452: * Compresses by zero-ing any value below a given percentile cut-off.
453: * @param p percentile cut-off, must be between 0 and 1.
454: * @param k the index of the coefficient array to compress.
455: **********************************************/
456: public void compress(double p, int k) {
457: coefs[k] = compress(coefs[k], p);
458: }
459:
460: /**************************************
461: * Compresses by zero-ing any value below a given percentile cut-off.
462: * @param v an array to compress.
463: * @param p percentile cut-off, must be between 0 and 1.
464: ***************************************/
465: public static double[] compress(double[] v, double p) {
466: if (p == 0)
467: return (v);
468: double[] ans = v;
469: double seuil = ArrayMath.percentile(ArrayMath.abs(ans), p);
470: for (int k = 0; k < ans.length; k++) {
471: if (Math.abs(ans[k]) <= seuil) {
472: ans[k] = 0;
473: }
474: }
475: return (ans);
476: }
477:
478: /*********************************************
479: * Denoises by zero-ing any value above a given cut-off.
480: * @param p cut-off.
481: **********************************************/
482: public void denoiseHard(double p) {
483: for (int k = 1; k < coefs.length; k++) {
484: coefs[k] = denoiseHard(coefs[k], p);
485: }
486: }
487:
488: /*********************************************
489: * Denoises by zero-ing any value above a given cut-off.
490: * @param p cut-off.
491: * @param k the index of the coefficient array to denoise.
492: **********************************************/
493: public void denoiseHard(double p, int k) {
494: coefs[k] = denoiseHard(coefs[k], p);
495: }
496:
497: /**************************************
498: * Denoises by zero-ing any value above a given cut-off.
499: * @param v an array to denoise.
500: * @param seuil cut-off/threshold.
501: ***************************************/
502: public static double[] denoiseHard(double[] v, double seuil) {
503: if (seuil < 0) {
504: throw new IllegalArgumentException(
505: "The cutoff value must be positive.");
506: }
507: double[] ans = v;
508: for (int k = 0; k < ans.length; k++) {
509: if (Math.abs(ans[k]) >= seuil) {
510: ans[k] = 0;
511: }
512: }
513: return (ans);
514: }
515:
516: /*********************************************
517: * Compresses by zero-ing any value below a given cut-off.
518: * @param p cut-off.
519: **********************************************/
520: public void compressHard(double p) {
521: for (int k = 1; k < coefs.length; k++) {
522: coefs[k] = compressHard(coefs[k], p);
523: }
524: }
525:
526: /*********************************************
527: * Compresses by zero-ing any value below a given cut-off.
528: * @param p cut-off.
529: * @param k the index of the coefficient array to compress.
530: **********************************************/
531: public void compressHard(double p, int k) {
532: coefs[k] = compressHard(coefs[k], p);
533: }
534:
535: /**************************************
536: * Compresses by zero-ing any value below a given cut-off.
537: * @param v an array to compress.
538: * @param seuil cut-off/threshold.
539: ***************************************/
540: public static double[] compressHard(double[] v, double seuil) {
541: if (seuil < 0) {
542: throw new IllegalArgumentException(
543: "The cutoff value must be positive.");
544: }
545: double[] ans = v;
546: for (int k = 0; k < ans.length; k++) {
547: if (Math.abs(ans[k]) <= seuil) {
548: ans[k] = 0;
549: }
550: }
551: return (ans);
552: }
553: }
|