VKEDCO: Vladimir Kulyukin's Education Coop: CDF22

Showing posts with label CDF22. Show all posts

Thursday, January 7, 2016

Computing Forward Ordered CDF(4, 4): Programmatic Note 7 on Ripples in Mathematics

Problem

This is my programmatic note 7 on "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo. One of the problems I had when I was working on CDF(2, 2) was the inverse transform. My values, when I was inversing the transformed signals, were off by 0.25 or some small amounts. I could not reconstruct the graphs in Fig. 4.10 and 4.11 in Ch. 4 in "Ripples in Mathematics." In general, my CDF implementation did not go nearly as smoothly as the Haar. I have documented my conceptual struggles and formalism modifications in my previous posts on CDF(2, 2) on this blog which you can find by searching it on "wavelets."

I started looking deeper into the coefficients for computing smoothed values and wavelets, i.e., H0, H1, G0, G1 for CDF(2, 2) and H0, H1, H2, H3, G0, G1, G2, G3, G4 for forward CDF(4, 4). What I found out was that these coefficients vary from publication to publication. Perhaps, I am still off on some of the finer mathematical points of Daubechies Wavelets and simply do not understand a formal abstraction (if one exists) where they are all the same. For those of you who are interested in digging deeper into these coefficients, here are the three references I looked at to better understand CDF(2, 2) and CDF(4, 4).

1) Y. Nievergelt. "Wavelets Made Easy";
2) G. Uytterhoeven, D. Roose, and A. Buttheel. "Wavelet Transforms Using the Lifting Scheme.";
3) G. Uytterhoeven, F. Van Wulpen, M. Jansen, D. Roose, and A. Bultheel. "WAILI: Wavelets with Integer Lifting."

These are great references and I keep going back to them and re-reading them, but they are very formal and lack programmatic details, especially for software developers working with imperative programming languages such as Java and C++. As luck would have it, I then found Ian Kaplan's wonderful site at www.bearcave.com. Ian has a whole bunch of great documents and C++ and Java programs on wavelets. After I read his explanation of CDF(4, 4), a lot of things became much clearer to me programmatically. Then I got back to implementing forward CDF(4, 4) and used Ian's coefficients. He credits Jensen and la Cour-Harbo in his implementation of CDF(4, 4). However, I could not find where in the book those coefficients are defined. I will keep looking for them though. Below is my implementation of forward CDF(4, 4). I plan to post my inverse implementation in a future blog entry. My Java source is here.

Forward CDF(4, 4) in Java

I started by defining constants for H0, H1, H2, H3 and G0, G1, G2, and G3.

public class CDF44 {
    static final double SQRT_OF_3 = Math.sqrt(3);
    static final double SQRT_OF_2 = Math.sqrt(2);
    static final double FOUR_SQRT_OF_2 = 4*SQRT_OF_2;

    // CDF(4,4) Forward signal coefficients
    static final double H0 = (1 + SQRT_OF_3)/FOUR_SQRT_OF_2;
    static final double H1 = (3 + SQRT_OF_3)/FOUR_SQRT_OF_2;
    static final double H2 = (3 - SQRT_OF_3)/FOUR_SQRT_OF_2;
    static final double H3 = (1 - SQRT_OF_3)/FOUR_SQRT_OF_2;

    // CDF(4, 4) Forward wavelet coefficients
    static final double G0 = H3;
    static final double G1 = -H2;
    static final double G2 = H1;
    static final double G3 = -H0;

}

Here is an implementation of forward CDF(4, 4) as a static method of CDF44. There is a dbg_flag variable that can be set to false when output is not needed. I generally prefer while-loops to for-loops. There are, and have always been, more intuitively obvious to me. Utils.isPowerOf2(N) is from my Utils class. It has a bunch of util methods, primarily numerical.

    public static void orderedDWT(double[] signal, boolean dbg_flag) {
        final int N = signal.length;
        if ( !Utils.isPowerOf2(N) ) return;
        if ( N < 4 ) return;
        int i, j, mid;
        double[] D4 = null;

        int numScalesToDo = Utils.powVal(N)-1;
        int currScale = 0;
        int signal_length = N;
        while ( signal_length >= 4 ) {
            mid = signal_length >> 1; // n / 2;
            if ( dbg_flag ) System.out.println("MID = " + mid);
            if ( dbg_flag ) System.out.println("signal_length   = " + signal_length);
            D4 = new double[signal_length]; // temporary array that saves the scalers and wavelets
            for(i = 0, j = 0; j < signal_length-3; i += 1, j += 2) {
                if ( dbg_flag ) {
                    final String cursig = "s^{" + (currScale+1) + "}_{" + (numScalesToDo-1) + "}";
                    final String prvsig = "s^{" + currScale + "}_{" + numScalesToDo + "}";
                    System.out.print("SCL: " + cursig + "[" + i + "]=" + "H0*" + prvsig + "[" + j + "]+H1*" + prvsig + "[" + (j+1) + "]+" +
                        "H2*" + prvsig + "[" + (j+2) + "]+" + "H3*" + prvsig + "[" + (j+3) + "]; " );
                    System.out.println("WVL: " + cursig + "[" + (mid+i) + "]=" + "G0*" + prvsig + "[" + j + "]+" + "G1*" + prvsig + "[" + (j+1) + "]+" +
                        "G2*" + prvsig + "[" + (j+2) + "]+" + "G3*" + prvsig + "[" + (j+3) + "]" );
                }
                // cdf44[i] is a scaled sample
                D4[i]     = H0*signal[j] + H1*signal[j+1] + H2*signal[j+2] + H3*signal[j+3];
                // cdf44[mid+i] is the corresponding wavelet for d4[i]
                D4[mid+i] = G0*signal[j] + G1*signal[j+1] + G2*signal[j+2] + G3*signal[j+3];
            }

            currScale     += 1;
            numScalesToDo -= 1;

            // cdf44[i] is a scaled sample with a mirror wrap-up
            D4[i]     = H0*signal[signal_length-2] + H1*signal[signal_length-1] + H2*signal[0] + H3*signal[1];
            // cdf44[mid+i] is the corresponding wavelet for d4[i]
            D4[mid+i] = G0*signal[signal_length-2] + G1*signal[signal_length-1] + G2*signal[0] + G3*signal[1];

            if ( dbg_flag ) {
                final String cursig = "s^{" + currScale + "}_{" + numScalesToDo + "}";
                final String prvsig = "s^{" + (currScale-1) + "}_{" + (numScalesToDo+1) + "}";
                System.out.print("SCL: " + cursig + "[" + i + "]=" + "H0*" + prvsig + "[" + (signal_length-2) + "]+H1*" + prvsig +
                                           "[" + (signal_length-1) + "]+" + "H2*" + prvsig + "[" + 0 + "]+" + "H3*" + prvsig + "[" + 1 + "]; " );
               System.out.println("WVL: " + cursig + "[" + (mid+i) + "]=" + "G0*" + prvsig + "[" + (signal_length-2) + "]+" +
                                               "G1*" + prvsig + "[" + (signal_length-1) + "]+" +
                                               "G2*" + prvsig + "[" + 0 + "]+" + "G3*" + prvsig + "[" + 1 + "]" );
            }

            System.arraycopy(D4, 0, signal, 0, D4.length);
            D4 = null;
            signal_length >>= 1; // signal_length gets halved at each iteration/scale
        }
    }

Tests

A few methods and constants to run some tests.

    public static void test_fwd_cdf44(double[] s, boolean dbg_flag) {
        double[] scopy = new double[s.length];
        System.arraycopy(s, 0, scopy, 0, s.length);
        System.out.print("Input: "); Utils.displaySample(scopy);
        CDF44.orderedDWT(s, dbg_flag);
        System.out.print("FWD CDF(4,4): "); Utils.displaySample(s);
        System.out.println();
    }

    static double[] a01a = {1, 2, 3, 4};
    static double[] a01b = {4, 3, 2, 1};
    static double[] a02a = {1, 2, 3, 4, 5, 6, 7, 8};
    static double[] a02b = {8, 7, 6, 5, 4, 3, 2, 1};

    static double[] a03a = {1, 1, 1, 1};
    static double[] a03b = {2, 2, 2, 2};
    static double[] a03c = {3, 3, 3, 3};
    static double[] a03d = {4, 4, 4, 4};

    static double[] a04a = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16};
    static double[] a04b = {16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1};

    public static void main(String[] args) {
        test_fwd_cdf44(a01a, true);
    }

Here is the output of main. You can plug in different arrays to test it. The value of the forward CDF(4, 4) is the last printed line.

Input: Sample: 4.0 3.0 2.0 1.0
MID           = 2
signal_length = 4
SCL: s^{1}_{0}[0]=H0*s^{0}_{1}[0]+H1*s^{0}_{1}[1]+H2*s^{0}_{1}[2]+H3*s^{0}_{1}[3]; WVL: s^{1}_{0}[2]=G0*s^{0}_{1}[0]+G1*s^{0}_{1}[1]+G2*s^{0}_{1}[2]+G3*s^{0}_{1}[3]
SCL: s^{1}_{0}[1]=H0*s^{0}_{1}[2]+H1*s^{0}_{1}[3]+H2*s^{0}_{1}[0]+H3*s^{0}_{1}[1]; WVL: s^{1}_{0}[3]=G0*s^{0}_{1}[2]+G1*s^{0}_{1}[3]+G2*s^{0}_{1}[0]+G3*s^{0}_{1}[1]
FWD CDF(4,4): Sample: 4.760278777324326 2.3107890345411484 -2.220446049250313E-16 1.4142135623730943

Tuesday, December 15, 2015

NCDF(S, 2, 2, L, 3): Example 3: S = [1, 1, 1, 1, 4, 4, 4, 4]

Vladimir Kulyukin

Problem

This post is a detailed example of computing individual elements of NCDF(S, 2, 2, L, 3) at each of the three scales on S = [1, 1, 1, 1, 4, 4, 4, 4]. NCDF(S, 2, 2, L, J) notation is defined here. This notation is based on the notation introduced in Chapters 3 and 4 in "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo.

NCDF(S, 2, 2, 1, 3): Scale L = 1

Figures 1 and 2 show the formulas for the elements of S and D at scale L = 1. Figure 3 shows an example of computing the first scale of the signal S = [1, 1, 1, 1, 4, 4, 4, 4].

Figure 1. Formulas for elements of S at scale L = 1

Figure 2. Formulas for elements of D at scale L = 1

Figure 3. Computing elements of S and D at scale L = 1

NCDF(S, 2, 2, 2, 3): Scale L = 2

Figure 4 shows the individual formulas for S and D at scale L = 2. Figure 5 shows the computation of individual elements of S and D of the sample S = [1, 1, 1, 1, 4, 4, 4, 4] at scale L = 2.


Figure 4. Formulas for elements of S and D at scale L = 2

Figure 5. Computing elements of S and D at scale L = 2

NCDF(S, 2, 2, 3, 3): Scale L = 3

Figure 6 shows the individual formulas for S and D at scale L = 3. Figure 7 shows the computation of individual elements of S and D of the sample S = [1, 1, 1, 1, 4, 4, 4, 4] at scale L = 3.

Figure 6. Formulas for elements of S and D at scale L = 3

Figure 7. Computing elements of S and D at scale L = 3

Summary of NCDF(S, 2, 2, 3, 3) at All Three Scales

Figure 8 summarizes the three scales of NCDF(S, 2, 2, 3, 3), where S = [1, 1, 1, 1, 4, 4, 4, 4].

Figure 8. Summary of three scales of NCDF(S, 2, 2, 3, 3) on S = [1, 1, 1, 1, 4, 4, 4, 4]

SCDF(S, 2, 2, L, 3): Example 3: S = [1, 1, 1, 1, 4, 4, 4, 4]

Vladimir Kulyukin

Problem

This post is a detailed example of computing individual elements of SCDF(S, 2, 2, L, 3) at each of the three scales on S = [1, 1, 1, 1, 4, 4, 4, 4]. SCDF(S, 2, 2, L, J) notation is defined here. This notation is based on the notation introduced in Chapters 3 and 4 in "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo.

SCDF(S, 2, 2, 1, 3): Scale L = 1

Figures 1 and 2 show the formulas for the elements of S and D at scale L = 1. Figure 3 shows an example of computing the first scale of the signal S = [1, 1, 1, 1, 4, 4, 4, 4].

Figure 1. Formulas for individual elements of S at scale L = 1

Figure 2. Formulas for individual elements of D at scale L = 1

Figure 3. Computing elements of S and D at scale L = 1

SCDF(S, 2, 2, 2, 3): Scale L = 2

Figure 4 shows the individual formulas for S and D at scale L = 2. Figure 5 shows the computation of individual elements of S and D of the sample S = [1, 1, 1, 1, 4, 4, 4, 4] at scale L = 2.

Figure 4. Formulas for individual elements of S and D at scale L = 2

Figure 5. Computing elements of S and D at scale L = 2

SCDF(S, 2, 2, 3, 3): Scale L = 3

Figure 6 shows the individual formulas for S and D at scale L = 3. Figure 7 shows the computation of individual elements of S and D of the sampleS = [1, 1, 1, 1, 4, 4, 4, 4] at scale L = 3.

Figure 6. Formulas for individual elements of S and D at scale L = 3

Figure 7. Computing elements of S and D at scale L = 3

Summary of SCDF(S, 2, 2, 3, 3) at All Three Scales

Figure 8 summarizes the three scales of SCDF(S, 2, 2, 3, 3), where S = S = [1, 1, 1, 1, 4, 4, 4, 4] .

Figure 8. Summary of SCDF(S, 2, 2, 3, 3) for three scales on S = [1, 1, 1, 1, 4, 4, 4, 4]

NCDF(S, 2, 2, L, 3): Example 2: S = [8, 7, 6, 5, 4, 3, 2, 1]

Vladimir Kulyukin

Problem

This post is a detailed example of computing individual elements of NCDF(S, 2, 2, L, 3) at each of the three scales on S = [8, 7, 6, 5, 4, 3, 2, 1]. NCDF(S, 2, 2, L, J) notation is defined here. This notation is based on the notation introduced in Chapters 3 and 4 in "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo.


Figure 1. Formulas for individual elements of S at scale L = 1

Figure 2. Formulas for individual elements of D at scale L = 1

Figure 3. Computing individual elements of S and D at scale L = 1

Figure 4. Formulas for individual elements of S and D at scale L = 2

Figure 5. Computing individual elements of S and D at scale L = 2


Figure 6. Formulas for individual elements of S and D at scale L = 3

Figure 7. Computing elements of S and D at scale L = 3

Summary of NCDF(S, 2, 2, 3, 3) at All Three Scales

Figure 8 summarizes the three scales of NCDF(S, 2, 2, 3, 3), where S = [8, 7, 6, 5, 4, 3, 2, 1].

Figure 8. Summary of NCDF(S, 2, 2, 3, 3) at all three scales on S = [8, 7, 6, 5, 4, 3, 2, 1]

SCDF(S, 2, 2, L, 3): Example 2: S = [8, 7, 6, 5, 4, 3, 2, 1]

Vladimir Kulyukin

Problem

This post is a detailed example of computing individual elements of SCDF(S, 2, 2, L, 3) at each of the three scales on S = [8, 7, 6, 5, 4, 3, 2, 1]. SCDF(S, 2, 2, L, J) notation is defined here. This notation is based on the notation introduced in Chapters 3 and 4 in "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo.

Figure 1. Formulas for elements of S at scale L = 1

Figure 2. Formulas for elements of D at scale L = 1

Figure 3. Computing elements of S and D at scale L = 1

Figure 4. Formulas for elements of S and D at scale L = 2

Figure 5. Computing elements of S and D at scale L = 2

Figure 6. Formulas for individual elements of S and D at scale L = 3

Figure 7. Computing elements of S and D at scale L = 3

Summary of SCDF(S, 2, 2, 3, 3) at All Three Scales

Figure 8 summarizes the three scales of SCDF(S, 2, 2, 3, 3), where S = [8, 7, 6, 5, 4, 3, 2, 1] .

Figure 8. Computing three scales of SCDF(S, 2, 2, 3, 3) for S = [8, 7, 6, 5, 4, 3, 2, 1]

NCDF(S, 2, 2, L, 3): Example 1: S = [1, 2, 3, 4, 5, 6, 7, 8]

Vladimir Kulyukin

Problem

This post is a detailed example of computing individual elements of NCDF(S, 2, 2, L, 3) at each of the three scales on S = [1, 2, 3, 4, 5, 6, 7, 8]. NCDF(S, 2, 2, L, J) notation is defined here. This notation is based on the notation introduced in Chapters 3 and 4 in "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo.

Figure 1. Formulas for individual elements of S at scale L = 1

Figure 2. Formulas for individual elements of D at scale L = 1

Figure 3. Computing elements of S and D at scale L = 1

Figure 4. Formulas for individual elements of S and D at scale L = 2

Figure 5. Computing elements of S and D at scale L = 2

Figure 6. Formulas for individual elements of S and D at scale L = 3

Figure 7. Computing elements of S and D at scale L = 3

Summary of NCDF(S, 2, 2, 3, 3) at All Three Scales

Figure 8 summarizes the three scales of NCDF(S, 2, 2, 3, 3), where S = [1, 2, 3, 4, 5, 6, 7, 8].

Figure 8. Summary of NCDF(S, 2, 2, 3, 3) at all three scales

SCDF(S, 2, 2, L, 3): Example 1: S = [1, 2, 3, 4, 5, 6, 7, 8]

Vladimir Kulyukin

Problem

This post is a detailed example of computing individual elements of SCDF(S, 2, 2, L, 3) at each of the three scales on S = [1, 2, 3, 4, 5, 6, 7, 8]. SCDF(S, 2, 2, L, J) notation is defined here. This notation is based on the notation introduced in Chapters 3 and 4 in "Ripples in Mathematics" by A. Jensen & A. la Cour-Harbo.

Figure 1. Formulas for individual elements of S at scale L = 1

Figure 2. Formulas for individual elements of D at scale L = 1

Figure 3. Computing individual elements of S and D at scale L = 1

Figure 4. Formulas for individual elements of S and D at scale L = 2

Figure 5. Computing individual elements of S and D of [1, 2, 3, 4, 5, 6, 7, 8] at scale L = 2

SCDF(S, 2, 2, 3, 3): Scale L = 3

Figure 6 shows the individual formulas for S and D at scale L = 3. Figure 8 shows the computation of individual elements of S and D of the sample [1, 2, 3, 4, 5, 6, 7, 8] at scale L = 3.

Figure 6. Formulas for individual elements of S and D at scale L = 3

Figure 7. Computing individual elements of S and D of [1, 2, 3, 4, 5, 6, 7, 8] at scale L = 3

Summary of SCDF(S, 2, 2, 3, 3) at All Three Scales

Figure 8 summarizes the three scales of SCDF(S, 2, 2, 3, 3), where S = [1, 2, 3, 4, 5, 6, 7, 8].

Figure 8. Summary of SCDF(S, 2, 2, 3, 3) Example 1