001/*******************************************************************************
002 * This software is provided as a supplement to the authors' textbooks on digital
003 * image processing published by Springer-Verlag in various languages and editions.
004 * Permission to use and distribute this software is granted under the BSD 2-Clause
005 * "Simplified" License (see http://opensource.org/licenses/BSD-2-Clause).
006 * Copyright (c) 2006-2023 Wilhelm Burger, Mark J. Burge. All rights reserved.
007 * Visit https://imagingbook.com for additional details.
008 ******************************************************************************/
009
010package imagingbook.spectral.dct;
011
012/**
013 * <p>
014 * Calculates the 1D discrete cosine transform (DFT) on {@code double} data.
015 * This is a naive implementation which calculates cosines repeatedly (i.e., does
016 * NOT use pre-calculated cosine tables), intended for demonstration purposes
017 * only. See Sec. 20.1 (Prog. 20.1) of [1] for additional details.
018 * </p>
019 * <p>
020 * Usage example:
021 * <pre>
022 * double[] data = {1, 2, 3, 4, 5, 3, 0};
023 * Dct1d.Double dct = new Dct1dSlow(data.length);
024 * dct.forward(data);
025 * dct.inverse(data);
026 * ... </pre>
027 * <p>
028 * [1] W. Burger, M.J. Burge, <em>Digital Image Processing &ndash; An Algorithmic
029 * Introduction</em>, 3rd ed, Springer (2022).
030 * </p>
031 * 
032 * @author WB
033 * @version 2022/10/23
034 * @see Dct1dDirect
035 * @see Dct1dFast
036 */
037public class Dct1dSlow implements Dct1d.Double {
038
039        private final double CM0 = 1.0 / Math.sqrt(2);
040        private final int M;
041        private final double[] tmp;
042
043        /**
044         * Constructor.
045         * @param M the data size
046         */
047        public Dct1dSlow(int M) {
048                this.M = M;
049                this.tmp = new double[M];
050        }
051        
052        @Override
053        public int getSize() {
054                return this.M;
055        }
056        
057        // ----------------------------------------------------------------------------
058
059        @Override
060        public void forward(double[] g) {
061                checkSize(g);
062                final double s = Math.sqrt(2.0 / M);
063                double[] G = tmp;
064                for (int m = 0; m < M; m++) {
065                        double cm = (m == 0) ? CM0 : 1.0;
066                        double sum = 0;
067                        for (int u = 0; u < M; u++) {
068                                double Phi = Math.PI * m * (2 * u + 1) / (2 * M);
069                                sum += g[u] * cm * Math.cos(Phi);
070                        }
071                        G[m] = s * sum;
072                }
073                System.arraycopy(G, 0, g, 0, M); // copy G -> g
074        }
075
076        @Override
077        public void inverse(double[] G) {
078                checkSize(G);
079                final double s = Math.sqrt(2.0 / M); //common scale factor
080                double[] g = tmp;
081                for (int u = 0; u < M; u++) {
082                        double sum = 0;
083                        for (int m = 0; m < M; m++) {
084                                double cm = (m == 0) ? CM0 : 1.0;
085                                double Phi = Math.PI * m * (2 * u + 1) / (2 * M);
086                                sum += G[m] * cm * Math.cos(Phi);
087                        }
088                        g[u] = s * sum;
089                }
090                System.arraycopy(g, 0, G, 0, M); // copy g -> G
091        }
092
093        // test example
094//      public static void main(String[] args) {
095//              PrintPrecision.set(6);
096//              double[] data = {1,2,3,4,5,3,0};
097//              System.out.println("Original data:      " + Matrix.toString(data));
098//
099//              Dct1d.Double dct = new Dct1dSlow(data.length);
100//              dct.forward(data);
101//              System.out.println("DCT spectrum:       " + Matrix.toString(data));
102//
103//              dct.inverse(data);
104//              System.out.println("Reconstructed data: " + Matrix.toString(data));
105//      }
106
107        //      Original data:      {1.000, 2.000, 3.000, 4.000, 5.000, 3.000, 0.000}
108        //      DCT of data:        {6.803, -0.361, -3.728, 1.692, -0.888, -0.083, 0.167}
109        //      Reconstructed data: {1.000, 2.000, 3.000, 4.000, 5.000, 3.000, -0.000}
110
111}