Merge branch 'feature/dct-iv' into 'master'

Add DCT/DST type IV

See merge request idf/esp-dsp!126

CHANGELOG.md

@@ -11,6 +11,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Bugfix for determinant calculation in mat.cpp 
 
 ### Added
+- Add DCT-IV and DST-IV 
 
 ### Removed
 

CMakeLists.txt

@@ -121,6 +121,8 @@ set(srcs            "modules/common/misc/dsps_pwroftwo.cpp"
                     "modules/fft/fixed/dsps_fft2r_sc16_arp4.S"
 
                     "modules/dct/float/dsps_dct_f32.c"
+                    "modules/dct/float/dsps_dctiv_f32.c"
+                    "modules/dct/float/dsps_dstiv_f32.c"
                     "modules/support/snr/float/dsps_snr_f32.cpp"
                     "modules/support/sfdr/float/dsps_sfdr_f32.cpp"
                     "modules/support/misc/dsps_d_gen.c"

modules/dct/float/dsps_dctiv_f32.c

@@ -0,0 +1,94 @@
+// Copyright 2018-2020 Espressif Systems (Shanghai) PTE LTD
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+#include "dsp_common.h"
+#include <math.h>
+
+#include "dsps_dct.h"
+#include "dsps_fft2r.h"
+
+esp_err_t dsps_dctiv_f32(float *data, int ndct)
+{
+    if (dsps_fft2r_initialized == 0) {
+        return ESP_ERR_DSP_REINITIALIZED;
+    }
+
+    float factor = M_PI / (ndct * 2);
+    float in1, in2, in3, in4;
+    for (int i = 0; i < ndct / 4; i++) {
+        in1 = data[i * 2 + 0];
+        in2 = data[i * 2 + 1];
+        in3 = data[ndct - i * 2 - 1];
+        in4 = data[ndct - i * 2 - 2];
+
+        data[i * 2 + 0] = (
+                              in1 * cos(factor * (i * 2 + 0))
+                              + in3 * cos(factor * ((ndct - i * 2)))
+                          );
+
+        data[i * 2 + 1] = (
+                              -in1 * sin(factor * (i * 2))
+                              + in3 * sin(factor * ((ndct - i * 2)))
+                          );
+
+        data[ndct - i * 2 + 0 - 2] = (
+                                         in2 * cos(factor * (i * 2 + 1 + 0.5) )
+                                         + in4 * cos(factor * ((ndct - i * 2 - 1) - 0.5) )
+                                     );
+
+        data[ndct - i * 2 + 1 - 2] = (
+                                         in2 * sin(factor * (i * 2 + 1))
+                                         + in4 * sin(-factor * ((ndct - i * 2 - 1)) )
+                                     );
+
+    }
+    esp_err_t error = ESP_OK;
+    error = dsps_fft2r_fc32(data, ndct / 2);
+    if (error != ESP_OK) {
+        return error;
+    }
+    error = dsps_bit_rev_fc32(data, ndct / 2);
+    if (error != ESP_OK) {
+        return error;
+    }
+
+    for (int i = 0; i < ndct / 4; i++) {
+        in1 = data[2 * i + 0];
+        in2 = data[2 * i + 1];
+
+        in3 = data[ndct - 2 * i - 2];
+        in4 = data[ndct - 2 * i - 1];
+
+        data[i * 2 + 0] = (
+                              in1 * cos(factor * (0 + i * 2))
+                              + in2 * sin(factor * (0 + i * 2))
+                          );
+
+        data[ndct  - i * 2 - 1] = (
+                                      in1 * cos(factor * (ndct - i * 2))
+                                      - in2 * sin(factor * (ndct - i * 2))
+                                  );
+
+        data[i * 2 + 1] = (
+                              in3 * cos(factor * (2 + i * 2))
+                              - in4 * sin(factor * (2 + i * 2))
+                          );
+
+        data[ndct  - i * 2 - 2] = (
+                                      in3 * cos(factor * (ndct - i * 2 - 2) )
+                                      + in4 * sin(factor * (ndct - i * 2 - 2) )
+                                  );
+    }
+    return ESP_OK;
+}

modules/dct/float/dsps_dstiv_f32.c

@@ -0,0 +1,98 @@
+// Copyright 2018-2020 Espressif Systems (Shanghai) PTE LTD
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//     http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+#include "dsp_common.h"
+#include <math.h>
+
+#include "dsps_dct.h"
+#include "dsps_fft2r.h"
+
+esp_err_t dsps_dstiv_f32(float *data, int ndst)
+{
+    if (dsps_fft2r_initialized == 0) {
+        return ESP_ERR_DSP_REINITIALIZED;
+    }
+
+    float in1, in2, in3, in4;
+    float factor = M_PI / (ndst);
+
+    for (int i = 0; i < ndst / 4; i++) {
+
+        in1 = data[2 * i + 0];
+        in2 = data[2 * i + 1];
+
+        in3 = data[ndst - 2 * i - 2];
+        in4 = data[ndst - 2 * i - 1];
+
+        data[i * 2 + 1] = (
+                              in1 * cos(factor * (i + 0))
+                              - in4 * sin(factor * ((ndst - i - 1)))
+                          );
+
+        data[i * 2 + 0] = (
+                              in1 * sin(factor * (i))
+                              - in4 * cos(factor * ((ndst - i - 1)))
+                          );
+
+        data[ndst - i * 2 - 2] = (
+                                     -in3 * cos(factor * (ndst - i - 1))
+                                     + in2 * sin(factor * (ndst - i - 1))
+                                 );
+
+        data[ndst - i * 2 - 1] = (
+                                     +in3 * sin(factor * (i + 1))
+                                     - in2 * cos(-factor * (i + 1))
+                                 );
+
+    }
+
+    esp_err_t error = ESP_OK;
+    error = dsps_fft2r_fc32(data, ndst / 2);
+    if (error != ESP_OK) {
+        return error;
+    }
+    error = dsps_bit_rev_fc32(data, ndst / 2);
+    if (error != ESP_OK) {
+        return error;
+    }
+
+    for (int i = 0; i < ndst / 4; i++) {
+        in1 = data[2 * i + 0];
+        in2 = data[2 * i + 1];
+
+        in3 = data[ndst - 2 * i - 2];
+        in4 = data[ndst - 2 * i - 1];
+
+        data[i * 2 + 0] = (
+                              in1 * cos(factor * (0 + i))
+                              + in2 * sin(factor * (0 + i))
+                          );
+
+        data[ndst - i * 2 - 2 + 1] = (
+                                         -in1 * cos(factor * (ndst / 2 - i))
+                                         + in2 * sin(factor * (ndst / 2 - i))
+                                     );
+
+        data[i * 2 + 1] = (
+                              -in3 * cos(factor * (1 + i))
+                              + in4 * sin(factor * (1 + i))
+                          );
+
+        data[ndst  - i * 2 - 2 + 0] = (
+                                          +in3 * cos(factor * (ndst / 2 - i - 1))
+                                          + in4 * sin(factor * (ndst / 2 - i - 1))
+                                      );
+    }
+    return ESP_OK;
+}

modules/dct/include/dsps_dct.h

@@ -26,7 +26,7 @@ extern "C"
 /**
  * @brief      DCT of radix 2, unscaled
  *
- * DCT type II of radix 2, unscaled
+ * Discrete Cosine Transform type II of radix 2, unscaled
  * Function is FFT based
  * The extension (_ansi) use ANSI C and could be compiled and run on any platform.
  * The extension (_ae32) is optimized for ESP32 chip.
@@ -41,14 +41,56 @@ extern "C"
  *      - One of the error codes from DSP library
  */
 esp_err_t dsps_dct_f32(float *data, int N);
+/**@}*/
+
+
+/**@{*/
+/**
+ * @brief      DCT of radix 2, type IV, unscaled
+ *
+ * Discrete Cosine Transform type IV of radix 2, unscaled
+ * Function is FFT based
+ * The extension (_ansi) use ANSI C and could be compiled and run on any platform.
+ * The extension (_ae32) is optimized for ESP32 chip.
+ *
+ * @param[inout] data: input/output array with size of N. An elements located: Re[0],Re[1], , ... Re[N-1]
+ *               result of DST will be stored to this array from 0...N-1.
+ *               Size of data array must be N
+ * @param[in] N: Size of DCT transform. Size of data array must be N
+ *
+ * @return
+ *      - ESP_OK on success
+ *      - One of the error codes from DSP library
+ */
+esp_err_t dsps_dctiv_f32(float *data, int N);
+/**@}*/
 
+/**@{*/
+/**
+ * @brief      DST of radix 2, type IV, unscaled
+ *
+ * Discrete Sine Transform type IV of radix 2, unscaled
+ * Function is FFT based
+ * The extension (_ansi) use ANSI C and could be compiled and run on any platform.
+ * The extension (_ae32) is optimized for ESP32 chip.
+ *
+ * @param[inout] data: input/output array with size of N*2. An elements located: Re[0],Re[1], , ... Re[N-1]
+ *               result of DST will be stored to this array from 0...N-1.
+ *               Size of data array must be N
+ * @param[in] N: Size of DST transform. Size of data array must be N
+ *
+ * @return
+ *      - ESP_OK on success
+ *      - One of the error codes from DSP library
+ */
+esp_err_t dsps_dstiv_f32(float *data, int N);
 /**@}*/
 
 /**@{*/
 /**
  * @brief      Inverce DCT of radix 2
  *
- * Inverce DCT type III of radix 2, unscaled
+ * Inverce Discrete Cosine Transform type II of radix 2, unscaled
  * Function is FFT based
  * The extension (_ansi) use ANSI C and could be compiled and run on any platform.
  * The extension (_ae32) is optimized for ESP32 chip.