Daubechies6


1   
2   package JSci.maths.wavelet.daubechies6;
3   
4   import JSci.maths.wavelet.*;
5   import JSci.maths.*;
6   import JSci.maths.wavelet.*;
7   
8   /******************************************
9   * Daubechies wavelets adapted to the
10  * interval by Meyer. Thanks to Pierre Vial
11  * for the filters.
12  * @author Daniel Lemire
13  *****************************************/
14  public final class Daubechies6 extends Multiresolution implements Filter, NumericalConstants {
15      protected final static int filtretype=10;
16      protected final static int minlength=20;
17      /****************************************
18    * This method is used to compute
19    * how the number of scaling functions
20    * changes from on scale to the other.
21    * Basically, if you have k scaling
22    * function and a Filter of type t, you'll
23    * have 2*k+t scaling functions at the
24    * next scale (dyadic case).
25    * Notice that this method assumes
26    * that one is working with the dyadic
27    * grid while the method "previousDimension"
28    * define in the interface "Filter" doesn't.
29      ******************************************/
30      public int getFilterType () {
31          return(filtretype);
32      }
33  
34      public  MultiscaleFunction primaryScaling(int n0, int k) {
35          return(new Scaling6(n0,k));
36      }
37  
38      public  MultiscaleFunction dualScaling(int n0, int k) {
39          return(new Scaling6(n0,k));
40      }
41      public  MultiscaleFunction primaryWavelet(int n0, int k) {
42          return(new Wavelet6(n0,k));
43      }
44      public  MultiscaleFunction dualWavelet(int n0, int k) {
45          return(new Wavelet6(n0,k));
46      }
47  
48      final static double[] vg={
49         0.0217847003266,
50         0.00493661237185,
51        - 0.166863215412,
52         - 0.0683231215866,
53         0.694457972958,
54         1.11389278393,
55         0.477904371333,
56         - 0.102724969862,
57         - 0.0297837512985,
58         0.06325056266,
59         0.00249992209279,
60         - 0.0110318675094};
61  
62      final static double[] v0temp={
63      0.8035811728E+00,
64     -0.5951951771E+00};
65      final static double[] v1temp={
66      0.2409444375E+00,
67      0.3253023901E+00,
68      0.8139624945E+00,
69     -0.4166403612E+00};
70      final static double[] v2temp={
71     -0.1713726513E+00,
72     -0.2313725672E+00,
73      0.2850349985E+00,
74      0.2770986033E+00,
75      0.5575844464E+00,
76     -0.6694543662E+00};
77      final static double[] v3temp={
78     -0.2876454491E-02,
79     -0.3883540664E-02,
80     -0.5655242446E-01,
81     -0.1151783354E+00,
82      0.5730520088E+00,
83      0.4076173301E+00,
84     -0.4430745754E+00,
85     -0.5409776766E+00};
86      final static double[] v4temp={
87      0.9214100804E-01,
88      0.1244008388E+00,
89      0.7306532844E-01,
90      0.2931573652E+00,
91     -0.1321263657E+00,
92     -0.2417688989E-01,
93     -0.3632781234E+00,
94      0.6792063579E-01,
95      0.1887648043E+00,
96     -0.8333645213E+00};
97      final static double[] v5temp={
98     -0.1061491551E+00,
99     -0.1433134304E+00,
100    -0.2622879500E+00,
101    -0.6856964195E+00,
102     0.2766563635E+00,
103    -0.8836692218E-01,
104     0.3328884755E+00,
105     0.1288342707E+00,
106     0.5008585553E-01,
107    -0.4619028205E+00,
108     0.8828610310E-02,
109    -0.3895963779E-01};
110     final static double[] v6temp={
111    -0.2023876977E-01,
112    -0.2732464068E-01,
113     0.4220601869E-01,
114     0.4941654301E-01,
115     0.3803171587E-01,
116     0.8472554296E-01,
117     0.1946516826E+00,
118    -0.6992838718E-01,
119     0.9438761468E+00,
120     0.1295249721E+00,
121    -0.7015973081E-01,
122     0.1735085659E+00,
123     0.4990121912E-02,
124    -0.2202083175E-01};
125     final static double[] v7temp={
126     0.5456011297E-01,
127     0.7366235691E-01,
128    -0.1344279366E+00,
129    -0.1735568822E+00,
130     0.3102530241E-02,
131    -0.1659152277E+00,
132    -0.3869255752E+00,
133     0.2453596290E+00,
134     0.2217987487E+00,
135     0.1874151629E+00,
136     0.7052410111E+00,
137    -0.3333463173E+00,
138    -0.4413043778E-01,
139     0.1047852370E+00,
140     0.3298337216E-02,
141    -0.1455518125E-01};
142     final static double[] v8temp={
143    -0.4678407876E-02,
144    -0.6316382647E-02,
145     0.3592640703E-03,
146    -0.6935341643E-02,
147     0.4493043023E-02,
148     0.4405183141E-02,
149     0.2315284183E-01,
150    -0.9374856771E-02,
151    -0.1151952189E+00,
152    -0.4165806865E-01,
153     0.4741171777E+00,
154     0.7976083337E+00,
155     0.3391664841E+00,
156    -0.7534748748E-01,
157    -0.2114679112E-01,
158     0.4508434214E-01,
159     0.1768528334E-02,
160    -0.7804311311E-02};
161     final static double[] v9temp={
162     0.4182788817E-04,
163     0.5647240556E-04,
164    -0.3442465609E-04,
165     0.1028288644E-05,
166     0.1389588842E-03,
167     0.7128147780E-04,
168    -0.1053053091E-03,
169     0.2899067329E-03,
170     0.1541349274E-01,
171     0.3547130473E-02,
172    -0.1181216165E+00,
173    -0.4823069885E-01,
174     0.4910659254E+00,
175     0.7876202163E+00,
176     0.3379287356E+00,
177    -0.7263467678E-01,
178    -0.2106028481E-01,
179     0.4472486618E-01,
180     0.1767711923E-02,
181    -0.7800708584E-02};
182 
183     final static double[] vd0temp={
184     0.6512594203E+00,
185     0.7588551689E+00};
186 
187     final static double[] vd1temp={
188    -0.6896282953E-01,
189     0.5918480130E-01,
190    -0.9848250307E+00,
191     0.1478544768E+00};
192 
193     final static double[] vd2temp={
194    -0.1667727153E-01,
195     0.1431265231E-01,
196     0.1090239805E+00,
197     0.7126760308E+00,
198    -0.1719730277E+00,
199     0.6709313982E+00};
200 
201     final static double[] vd3temp={
202     0.9785571254E-01,
203    -0.8398105097E-01,
204     0.6834310191E-02,
205     0.1247808574E+00,
206     0.7236789907E-01,
207    -0.1108818957E+00,
208     0.1804007431E+00,
209     0.9579551851E+00};
210 
211     final static double[] vd4temp={
212    -0.5378738739E+00,
213     0.4616103858E+00,
214     0.2733949616E-01,
215    -0.2535533852E+00,
216     0.3484637175E+00,
217     0.3309876136E+00,
218    -0.2950887461E+00,
219     0.1958019515E+00,
220     0.6099724771E-01,
221     0.2691738861E+00};
222 
223     final static double[] vd5temp={
224    -0.2335808916E+00,
225     0.2004621729E+00,
226     0.5067312783E-01,
227     0.1483313327E+00,
228     0.4193562472E+00,
229    -0.6838768540E-01,
230     0.8092621294E+00,
231    -0.1702431983E+00,
232    -0.5881963586E-01,
233    -0.1104006534E+00,
234     0.5469174586E-02,
235     0.2413483588E-01};
236 
237     final static double[] vd6temp={
238     0.7905558067E-01,
239    -0.6784653218E-01,
240    -0.2413743535E-02,
241     0.4795434784E-01,
242    -0.2979345169E-01,
243    -0.5477003063E-01,
244     0.1575224560E+00,
245    -0.5400590327E-01,
246     0.8600210456E+00,
247     0.4427411664E+00,
248    -0.5350546863E-01,
249    -0.1279990455E+00,
250     0.3964078767E-02,
251     0.1749302183E-01};
252 
253     final static double[] vd7temp={
254     0.4472356538E-01,
255    -0.3838234812E-01,
256    -0.1245771188E-02,
257     0.2792643621E-01,
258    -0.2285669131E-01,
259    -0.3338969826E-01,
260     0.4201463682E-01,
261    -0.2161240035E-01,
262    -0.6527800067E-01,
263     0.3288435947E+00,
264     0.7890314114E+00,
265     0.4902682296E+00,
266    -0.4839771981E-01,
267    -0.1181080824E+00,
268     0.3500297801E-02,
269     0.1544640999E-01};
270 
271     final static double[] vd8temp={
272    -0.2366701878E-04,
273     0.2031134472E-04,
274    -0.3373738591E-03,
275    -0.2266339823E-02,
276    -0.2165041937E-02,
277     0.1906206337E-02,
278    -0.7087323880E-02,
279     0.2020684368E-02,
280     0.4460763770E-01,
281    -0.2107187698E-01,
282    -0.7263524048E-01,
283     0.3379418668E+00,
284     0.7876410397E+00,
285     0.4910540301E+00,
286    -0.4831177914E-01,
287    -0.1179902457E+00,
288     0.3490713064E-02,
289     0.1540411365E-01};
290 
291     final static double[] vd9temp={
292    -0.1333785556E-05,
293     0.1144672190E-05,
294     0.4933806493E-06,
295     0.2205985224E-05,
296    -0.3259284258E-06,
297    -0.2564526296E-05,
298    -0.1197305244E-05,
299    -0.1010124086E-06,
300    -0.7800395608E-02,
301     0.1767376409E-02,
302     0.4472504864E-01,
303    -0.2106027212E-01,
304    -0.7263753168E-01,
305     0.3379294037E+00,
306     0.7876411418E+00,
307     0.4910559454E+00,
308    -0.4831174259E-01,
309    -0.1179901111E+00,
310     0.3490712084E-02,
311     0.1540410933E-01};
312 
313     final static double[] v0=ArrayMath.scalarMultiply(SQRT2,v0temp);
314     final static double[] v1=ArrayMath.scalarMultiply(SQRT2,v1temp);
315     final static double[] v2=ArrayMath.scalarMultiply(SQRT2,v2temp);
316     final static double[] v3=ArrayMath.scalarMultiply(SQRT2,v3temp);
317     final static double[] v4=ArrayMath.scalarMultiply(SQRT2,v4temp);
318     final static double[] v5=ArrayMath.scalarMultiply(SQRT2,v5temp);
319     final static double[] v6=ArrayMath.scalarMultiply(SQRT2,v6temp);
320     final static double[] v7=ArrayMath.scalarMultiply(SQRT2,v7temp);
321     final static double[] v8=ArrayMath.scalarMultiply(SQRT2,v8temp);
322     final static double[] v9=ArrayMath.scalarMultiply(SQRT2,v9temp);
323 
324     final static double[] vd0=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd0temp));
325     final static double[] vd1=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd1temp));
326     final static double[] vd2=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd2temp));
327     final static double[] vd3=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd3temp));
328     final static double[] vd4=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd4temp));
329     final static double[] vd5=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd5temp));
330     final static double[] vd6=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd6temp));
331     final static double[] vd7=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd7temp));
332     final static double[] vd8=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd8temp));
333     final static double[] vd9=ArrayMath.invert(ArrayMath.scalarMultiply(SQRT2,vd9temp));
334 
335 
336     /********************************************
337     * On d�finit ici le filtre comme tel par le
338     * vecteur phvg (filtre passe-haut).
339     *********************************************/
340     final static double[] vgtemp=ArrayMath.scalarMultiply(1.0/SQRT2,vg);
341     final static double[] phvg=WaveletMath.lowToHigh(vgtemp);
342     final static double[] phv0={
343    -0.4056098474E+00,
344    -0.5476194186E+00,
345     0.3345888248E+00,
346    -0.8469145514E-02,
347     0.2860429450E-02,
348     0.4344317225E+00,
349    -0.1793040681E+00,
350     0.4501354334E+00,
351    -0.7823562605E-02,
352    -0.2169451337E-03,
353     0.3443284776E-06,
354    -0.1518815163E-05};
355     final static double[] phv1={
356    -0.2714331420E+00,
357    -0.3664656082E+00,
358     0.8617650411E-01,
359    -0.2747405494E+00,
360    -0.5148935461E+00,
361    -0.3097406490E+00,
362    -0.2206953036E+00,
363    -0.5444907507E+00,
364     0.4953215623E-01,
365    -0.2013924147E-01,
366    -0.1989883688E-02,
367     0.4859370755E-02,
368     0.1437935110E-03,
369    -0.6345441585E-03};
370     final static double[] phv2={
371    -0.7388896229E-01,
372    -0.9975850151E-01,
373     0.1991481291E+00,
374     0.2684432545E+00,
375    -0.1045213665E-01,
376     0.2405920775E+00,
377     0.5401581689E+00,
378    -0.3490563093E+00,
379    -0.6165711214E-01,
380    -0.1943705670E+00,
381     0.5024780923E+00,
382    -0.3169028187E+00,
383    -0.3084242027E-01,
384     0.7560970841E-01,
385     0.2218060423E-02,
386    -0.9788044508E-02};
387     final static double[] phv3={
388     0.2329741271E-02,
389     0.3145415646E-02,
390    -0.1484555312E-03,
391     0.3513131058E-02,
392    -0.2412185766E-02,
393    -0.2301643677E-02,
394    -0.1162880173E-01,
395     0.4467413223E-02,
396     0.4333253688E-01,
397     0.1764767584E-01,
398    -0.6404967457E-01,
399    -0.3432838094E+00,
400     0.7870194497E+00,
401    -0.4897071000E+00,
402    -0.4826799594E-01,
403     0.1178085687E+00,
404     0.3490290194E-02,
405    -0.1540224757E-01};
406     final static double[] phv4={
407    -0.2118182730E-04,
408    -0.2859787559E-04,
409     0.1743279790E-04,
410    -0.5207298811E-06,
411    -0.7036938958E-04,
412    -0.3609725256E-04,
413     0.5332707089E-04,
414    -0.1468100425E-03,
415    -0.7805460453E-02,
416    -0.1796282458E-02,
417     0.4479149901E-01,
418     0.2101925260E-01,
419    -0.7264258737E-01,
420    -0.3379188415E+00,
421     0.7876414895E+00,
422    -0.4910573861E+00,
423    -0.4831174654E-01,
424     0.1179901294E+00,
425     0.3490712054E-02,
426    -0.1540410920E-01};
427     final static double[] phvd0temp={
428     0.4378023718E+00,
429    -0.3757277153E+00,
430    -0.6873475431E-01,
431    -0.1032247812E+00,
432     0.7266439214E+00,
433     0.3259673475E+00,
434    -0.1186803367E+00,
435    -0.5853826628E-01,
436    -0.3032048461E-03,
437    -0.4670847828E-04,
438     0.1703850068E-05,
439     0.7518930104E-05};
440     final static double[] phvd1temp={
441    -0.9713836045E-01,
442     0.8336541004E-01,
443     0.1006861439E+00,
444     0.5919694503E+00,
445     0.3702101657E+00,
446    -0.5544630394E+00,
447    -0.4149317683E+00,
448    -0.7460178119E-01,
449     0.2790428530E-01,
450     0.5912930318E-02,
451    -0.1712035030E-02,
452    -0.3636645827E-02,
453     0.1436692465E-03,
454     0.6339957943E-03};
455     final static double[] phvd2temp={
456    -0.8958975532E-01,
457     0.7688709850E-01,
458     0.1185119667E-02,
459    -0.6467017684E-01,
460     0.3743730652E-01,
461     0.7423011440E-01,
462    -0.1113332259E+00,
463     0.5103740064E-01,
464     0.4816519926E+00,
465    -0.7799800339E+00,
466     0.3345221160E+00,
467     0.7239621136E-01,
468    -0.2085354320E-01,
469    -0.4433542083E-01,
470     0.1748540164E-02,
471     0.7716105825E-02};
472     final static double[] phvd3temp={
473     0.4166215409E-04,
474    -0.3575500495E-04,
475     0.6680919371E-03,
476     0.4483753011E-02,
477     0.4274088494E-02,
478    -0.3773957888E-02,
479     0.1399084946E-01,
480    -0.3990648231E-02,
481    -0.1177575130E+00,
482     0.4833340703E-01,
483     0.4910525706E+00,
484    -0.7876666869E+00,
485     0.3379295212E+00,
486     0.7264103605E-01,
487    -0.2106021362E-01,
488    -0.4472460633E-01,
489     0.1767709932E-02,
490     0.7800699800E-02};
491     final static double[] phvd4temp={
492     0.2633834884E-05,
493    -0.2260391508E-05,
494    -0.9742819475E-06,
495    -0.4356173277E-05,
496     0.6436128145E-06,
497     0.5064186603E-05,
498     0.2364326456E-05,
499     0.1994698649E-06,
500     0.1540349180E-01,
501    -0.3490049659E-02,
502    -0.1179904012E+00,
503     0.4831170232E-01,
504     0.4910559595E+00,
505    -0.7876411054E+00,
506     0.3379294202E+00,
507     0.7263751587E-01,
508    -0.2106029251E-01,
509    -0.4472490177E-01,
510     0.1767711864E-02,
511     0.7800708325E-02};
512 
513     final static double[] phvd0=ArrayMath.invert(phvd0temp);
514     final static double[] phvd1=ArrayMath.invert(phvd1temp);
515     final static double[] phvd2=ArrayMath.invert(phvd2temp);
516     final static double[] phvd3=ArrayMath.invert(phvd3temp);
517     final static double[] phvd4=ArrayMath.invert(phvd4temp);
518 
519     /****************************************
520   * This method return the number of "scaling"
521   * functions at the previous scale given a
522   * number of scaling functions. The answer
523   * is always smaller than the provided value
524   * (about half since this is a dyadic
525   * implementation). This relates to the same idea
526   * as the "Filter type". It is used by
527   * the interface "Filter".
528     *****************************************/
529     public int previousDimension (int k) {
530         return(Cascades.previousDimension(filtretype,k));
531 
532     }
533 
534 
535     public Daubechies6 () {}
536     /****************************************
537     * This is the implementation of the lowpass
538   * Filter. It is used by the interface
539   * "Filter". Lowpass filters are normalized
540   * so that they preserve constants away from
541   * the boundaries.
542     *****************************************/
543     public double[] lowpass (double[] v, double[] param) {
544         return(lowpass(v));
545     }
546     /****************************************
547     * This is the implementation of the highpass
548   * Filter. It is used by the interface
549   * "Filter". Highpass filters are normalized
550   * in order to get L2 orthonormality of the
551   * resulting wavelets (when it applies).
552   * See the class DiscreteHilbertSpace for
553   * an implementation of the L2 integration.
554     *****************************************/
555     public double[] highpass (double[] v, double[] param) {
556         return(highpass(v));
557     }
558     /****************************************
559     * This is the implementation of the lowpass
560   * Filter. It is used by the interface
561   * "Filter". Lowpass filters are normalized
562   * so that they preserve constants away from
563   * the boundaries.
564     *****************************************/
565     public double[] lowpass (double[] gete) {
566         if(gete.length<minlength) {
567             throw new IllegalScalingException("The array is not long enough : "+gete.length+" < "+minlength);
568         }
569         double[] sortie=new double[2*gete.length-filtretype];
570         int dl0=gete.length-1;
571         for(int k=10;k<=dl0-10;k++) {
572             for(int L=-6;L<6;L++){
573                 sortie[2*k+L-4]+=vg[L+6]*gete[k];
574             }
575         }
576         sortie=ArrayMath.add(sortie,gete[0],v0,0);
577         sortie=ArrayMath.add(sortie,gete[1],v1,0);
578         sortie=ArrayMath.add(sortie,gete[2],v2,0);
579         sortie=ArrayMath.add(sortie,gete[3],v3,0);
580         sortie=ArrayMath.add(sortie,gete[4],v4,0);
581         sortie=ArrayMath.add(sortie,gete[5],v5,0);
582         sortie=ArrayMath.add(sortie,gete[6],v6,0);
583         sortie=ArrayMath.add(sortie,gete[7],v7,0);
584         sortie=ArrayMath.add(sortie,gete[8],v8,0);
585         sortie=ArrayMath.add(sortie,gete[9],v9,0);
586         int p0=sortie.length-vd0.length;
587         int p1=sortie.length-vd1.length;
588         int p2=sortie.length-vd2.length;
589         int p3=sortie.length-vd3.length;
590         int p4=sortie.length-vd4.length;
591         int p5=sortie.length-vd5.length;
592         int p6=sortie.length-vd6.length;
593         int p7=sortie.length-vd7.length;
594         int p8=sortie.length-vd8.length;
595         int p9=sortie.length-vd9.length;
596         sortie=ArrayMath.add(sortie,gete[dl0],vd0,p0);
597         sortie=ArrayMath.add(sortie,gete[dl0-1],vd1,p1);
598         sortie=ArrayMath.add(sortie,gete[dl0-2],vd2,p2);
599         sortie=ArrayMath.add(sortie,gete[dl0-3],vd3,p3);
600         sortie=ArrayMath.add(sortie,gete[dl0-4],vd4,p4);
601         sortie=ArrayMath.add(sortie,gete[dl0-5],vd5,p5);
602         sortie=ArrayMath.add(sortie,gete[dl0-6],vd6,p6);
603         sortie=ArrayMath.add(sortie,gete[dl0-7],vd7,p7);
604         sortie=ArrayMath.add(sortie,gete[dl0-8],vd8,p8);
605         sortie=ArrayMath.add(sortie,gete[dl0-9],vd9,p9);
606         return(sortie);
607     }
608 
609     /****************************************
610     * This is the implementation of the highpass
611   * Filter. It is used by the interface
612   * "Filter". Highpass filters are normalized
613   * in order to get L2 orthonormality of the
614   * resulting wavelets (when it applies).
615   * See the class DiscreteHilbertSpace for
616   * an implementation of the L2 integration.
617     *****************************************/
618     public double[] highpass(double[] gete) {
619         if(gete.length<10) {
620             throw new IllegalScalingException("The array is not long enough : "+gete.length+" < 10");
621         }
622         double[] sortie=new double[2*gete.length+filtretype];
623         int dl0=gete.length-1;
624         for(int k=5;k<=dl0-5;k++) {
625             for(int L=-6;L<6;L++){
626                 sortie[2*k+L + 6 ]+=phvg[L+6]*gete[k];
627             }
628         }
629         sortie=ArrayMath.add(sortie,gete[0],phv0,0);
630         int p0=sortie.length-phvd0.length;
631         sortie=ArrayMath.add(sortie,gete[dl0],phvd0,p0);
632         sortie=ArrayMath.add(sortie,gete[1],phv1,0);
633         int p1=sortie.length-phvd1.length;
634         sortie=ArrayMath.add(sortie,gete[dl0-1],phvd1,p1);
635         sortie=ArrayMath.add(sortie,gete[2],phv2,0);
636         int p2=sortie.length-phvd2.length;
637         sortie=ArrayMath.add(sortie,gete[dl0-2],phvd2,p2);
638         sortie=ArrayMath.add(sortie,gete[3],phv3,0);
639         int p3=sortie.length-phvd3.length;
640         sortie=ArrayMath.add(sortie,gete[dl0-3],phvd3,p3);
641         sortie=ArrayMath.add(sortie,gete[4],phv4,0);
642         int p4=sortie.length-phvd4.length;
643         sortie=ArrayMath.add(sortie,gete[dl0-4],phvd4,p4);
644 
645         return(sortie);
646     }
647 
648 
649 
650     public double[] evalScaling (int n0, int k, int j1) {
651         return(Cascades.evalScaling(this,n0,j1,k));
652     }
653 
654 
655     public double[] evalWavelet (int n0, int k, int j1) {
656         return(Cascades.evalWavelet(this,filtretype,n0,j1,k));
657     }
658 
659 }
660
A to Z: JavaDoc & Examples Daily Java News & Articles Open Source Projects Open Source Codes Free Computer Books Remove Frame
Free Books Free Magazines
Popular Tags