1 /* ----------------------------------------------------------------------
2 * Project: CMSIS DSP Library
3 * Title: arm_dct4_q15.c
4 * Description: Processing function of DCT4 & IDCT4 Q15
5 *
6 * $Date: 23 April 2021
7 * $Revision: V1.9.0
8 *
9 * Target Processor: Cortex-M and Cortex-A cores
10 * -------------------------------------------------------------------- */
11 /*
12 * Copyright (C) 2010-2021 ARM Limited or its affiliates. All rights reserved.
13 *
14 * SPDX-License-Identifier: Apache-2.0
15 *
16 * Licensed under the Apache License, Version 2.0 (the License); you may
17 * not use this file except in compliance with the License.
18 * You may obtain a copy of the License at
19 *
20 * www.apache.org/licenses/LICENSE-2.0
21 *
22 * Unless required by applicable law or agreed to in writing, software
23 * distributed under the License is distributed on an AS IS BASIS, WITHOUT
24 * WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
25 * See the License for the specific language governing permissions and
26 * limitations under the License.
27 */
28
29 #include "dsp/transform_functions.h"
30
31 /**
32 @addtogroup DCT4_IDCT4
33 @{
34 */
35
36 /**
37 @brief Processing function for the Q15 DCT4/IDCT4.
38 @param[in] S points to an instance of the Q15 DCT4 structure.
39 @param[in] pState points to state buffer.
40 @param[in,out] pInlineBuffer points to the in-place input and output buffer.
41 @return none
42
43 @par Input an output formats
44 Internally inputs are downscaled in the RFFT process function to avoid overflows.
45 Number of bits downscaled, depends on the size of the transform. The input and output
46 formats for different DCT sizes and number of bits to upscale are mentioned in the table below:
47
48 \image html dct4FormatsQ15Table.gif
49 */
50
arm_dct4_q15(const arm_dct4_instance_q15 * S,q15_t * pState,q15_t * pInlineBuffer)51 void arm_dct4_q15(
52 const arm_dct4_instance_q15 * S,
53 q15_t * pState,
54 q15_t * pInlineBuffer)
55 {
56 const q15_t *weights = S->pTwiddle; /* Pointer to the Weights table */
57 const q15_t *cosFact = S->pCosFactor; /* Pointer to the cos factors table */
58 q15_t *pS1, *pS2, *pbuff; /* Temporary pointers for input buffer and pState buffer */
59 q15_t in; /* Temporary variable */
60 uint32_t i; /* Loop counter */
61
62
63 /* DCT4 computation involves DCT2 (which is calculated using RFFT)
64 * along with some pre-processing and post-processing.
65 * Computational procedure is explained as follows:
66 * (a) Pre-processing involves multiplying input with cos factor,
67 * r(n) = 2 * u(n) * cos(pi*(2*n+1)/(4*n))
68 * where,
69 * r(n) -- output of preprocessing
70 * u(n) -- input to preprocessing(actual Source buffer)
71 * (b) Calculation of DCT2 using FFT is divided into three steps:
72 * Step1: Re-ordering of even and odd elements of input.
73 * Step2: Calculating FFT of the re-ordered input.
74 * Step3: Taking the real part of the product of FFT output and weights.
75 * (c) Post-processing - DCT4 can be obtained from DCT2 output using the following equation:
76 * Y4(k) = Y2(k) - Y4(k-1) and Y4(-1) = Y4(0)
77 * where,
78 * Y4 -- DCT4 output, Y2 -- DCT2 output
79 * (d) Multiplying the output with the normalizing factor sqrt(2/N).
80 */
81
82 /*-------- Pre-processing ------------*/
83 /* Multiplying input with cos factor i.e. r(n) = 2 * x(n) * cos(pi*(2*n+1)/(4*n)) */
84 arm_mult_q15 (pInlineBuffer, cosFact, pInlineBuffer, S->N);
85 arm_shift_q15 (pInlineBuffer, 1, pInlineBuffer, S->N);
86
87 /* ----------------------------------------------------------------
88 * Step1: Re-ordering of even and odd elements as
89 * pState[i] = pInlineBuffer[2*i] and
90 * pState[N-i-1] = pInlineBuffer[2*i+1] where i = 0 to N/2
91 ---------------------------------------------------------------------*/
92
93 /* pS1 initialized to pState */
94 pS1 = pState;
95
96 /* pS2 initialized to pState+N-1, so that it points to the end of the state buffer */
97 pS2 = pState + (S->N - 1U);
98
99 /* pbuff initialized to input buffer */
100 pbuff = pInlineBuffer;
101
102
103 #if defined (ARM_MATH_LOOPUNROLL)
104
105 /* Initializing the loop counter to N/2 >> 2 for loop unrolling by 4 */
106 i = S->Nby2 >> 2U;
107
108 /* First part of the processing with loop unrolling. Compute 4 outputs at a time.
109 ** a second loop below computes the remaining 1 to 3 samples. */
110 do
111 {
112 /* Re-ordering of even and odd elements */
113 /* pState[i] = pInlineBuffer[2*i] */
114 *pS1++ = *pbuff++;
115 /* pState[N-i-1] = pInlineBuffer[2*i+1] */
116 *pS2-- = *pbuff++;
117
118 *pS1++ = *pbuff++;
119 *pS2-- = *pbuff++;
120
121 *pS1++ = *pbuff++;
122 *pS2-- = *pbuff++;
123
124 *pS1++ = *pbuff++;
125 *pS2-- = *pbuff++;
126
127 /* Decrement loop counter */
128 i--;
129 } while (i > 0U);
130
131 /* pbuff initialized to input buffer */
132 pbuff = pInlineBuffer;
133
134 /* pS1 initialized to pState */
135 pS1 = pState;
136
137 /* Initializing the loop counter to N/4 instead of N for loop unrolling */
138 i = S->N >> 2U;
139
140 /* Processing with loop unrolling 4 times as N is always multiple of 4.
141 * Compute 4 outputs at a time */
142 do
143 {
144 /* Writing the re-ordered output back to inplace input buffer */
145 *pbuff++ = *pS1++;
146 *pbuff++ = *pS1++;
147 *pbuff++ = *pS1++;
148 *pbuff++ = *pS1++;
149
150 /* Decrement the loop counter */
151 i--;
152 } while (i > 0U);
153
154
155 /* ---------------------------------------------------------
156 * Step2: Calculate RFFT for N-point input
157 * ---------------------------------------------------------- */
158 /* pInlineBuffer is real input of length N , pState is the complex output of length 2N */
159 arm_rfft_q15 (S->pRfft, pInlineBuffer, pState);
160
161 /*----------------------------------------------------------------------
162 * Step3: Multiply the FFT output with the weights.
163 *----------------------------------------------------------------------*/
164 arm_cmplx_mult_cmplx_q15 (pState, weights, pState, S->N);
165
166 /* The output of complex multiplication is in 3.13 format.
167 * Hence changing the format of N (i.e. 2*N elements) complex numbers to 1.15 format by shifting left by 2 bits. */
168 arm_shift_q15 (pState, 2, pState, S->N * 2);
169
170 /* ----------- Post-processing ---------- */
171 /* DCT-IV can be obtained from DCT-II by the equation,
172 * Y4(k) = Y2(k) - Y4(k-1) and Y4(-1) = Y4(0)
173 * Hence, Y4(0) = Y2(0)/2 */
174 /* Getting only real part from the output and Converting to DCT-IV */
175
176 /* Initializing the loop counter to N >> 2 for loop unrolling by 4 */
177 i = (S->N - 1U) >> 2U;
178
179 /* pbuff initialized to input buffer. */
180 pbuff = pInlineBuffer;
181
182 /* pS1 initialized to pState */
183 pS1 = pState;
184
185 /* Calculating Y4(0) from Y2(0) using Y4(0) = Y2(0)/2 */
186 in = *pS1++ >> 1U;
187 /* input buffer acts as inplace, so output values are stored in the input itself. */
188 *pbuff++ = in;
189
190 /* pState pointer is incremented twice as the real values are located alternatively in the array */
191 pS1++;
192
193 /* First part of the processing with loop unrolling. Compute 4 outputs at a time.
194 ** a second loop below computes the remaining 1 to 3 samples. */
195 do
196 {
197 /* Calculating Y4(1) to Y4(N-1) from Y2 using equation Y4(k) = Y2(k) - Y4(k-1) */
198 /* pState pointer (pS1) is incremented twice as the real values are located alternatively in the array */
199 in = *pS1++ - in;
200 *pbuff++ = in;
201 /* points to the next real value */
202 pS1++;
203
204 in = *pS1++ - in;
205 *pbuff++ = in;
206 pS1++;
207
208 in = *pS1++ - in;
209 *pbuff++ = in;
210 pS1++;
211
212 in = *pS1++ - in;
213 *pbuff++ = in;
214 pS1++;
215
216 /* Decrement the loop counter */
217 i--;
218 } while (i > 0U);
219
220 /* If the blockSize is not a multiple of 4, compute any remaining output samples here.
221 ** No loop unrolling is used. */
222 i = (S->N - 1U) % 0x4U;
223
224 while (i > 0U)
225 {
226 /* Calculating Y4(1) to Y4(N-1) from Y2 using equation Y4(k) = Y2(k) - Y4(k-1) */
227 /* pState pointer (pS1) is incremented twice as the real values are located alternatively in the array */
228 in = *pS1++ - in;
229 *pbuff++ = in;
230
231 /* points to the next real value */
232 pS1++;
233
234 /* Decrement loop counter */
235 i--;
236 }
237
238
239 /*------------ Normalizing the output by multiplying with the normalizing factor ----------*/
240
241 /* Initializing the loop counter to N/4 instead of N for loop unrolling */
242 i = S->N >> 2U;
243
244 /* pbuff initialized to the pInlineBuffer(now contains the output values) */
245 pbuff = pInlineBuffer;
246
247 /* Processing with loop unrolling 4 times as N is always multiple of 4. Compute 4 outputs at a time */
248 do
249 {
250 /* Multiplying pInlineBuffer with the normalizing factor sqrt(2/N) */
251 in = *pbuff;
252 *pbuff++ = ((q15_t) (((q31_t) in * S->normalize) >> 15));
253
254 in = *pbuff;
255 *pbuff++ = ((q15_t) (((q31_t) in * S->normalize) >> 15));
256
257 in = *pbuff;
258 *pbuff++ = ((q15_t) (((q31_t) in * S->normalize) >> 15));
259
260 in = *pbuff;
261 *pbuff++ = ((q15_t) (((q31_t) in * S->normalize) >> 15));
262
263 /* Decrement loop counter */
264 i--;
265 } while (i > 0U);
266
267
268 #else
269
270 /* Initializing the loop counter to N/2 */
271 i = S->Nby2;
272
273 do
274 {
275 /* Re-ordering of even and odd elements */
276 /* pState[i] = pInlineBuffer[2*i] */
277 *pS1++ = *pbuff++;
278 /* pState[N-i-1] = pInlineBuffer[2*i+1] */
279 *pS2-- = *pbuff++;
280
281 /* Decrement the loop counter */
282 i--;
283 } while (i > 0U);
284
285 /* pbuff initialized to input buffer */
286 pbuff = pInlineBuffer;
287
288 /* pS1 initialized to pState */
289 pS1 = pState;
290
291 /* Initializing the loop counter */
292 i = S->N;
293
294 do
295 {
296 /* Writing the re-ordered output back to inplace input buffer */
297 *pbuff++ = *pS1++;
298
299 /* Decrement the loop counter */
300 i--;
301 } while (i > 0U);
302
303
304 /* ---------------------------------------------------------
305 * Step2: Calculate RFFT for N-point input
306 * ---------------------------------------------------------- */
307 /* pInlineBuffer is real input of length N , pState is the complex output of length 2N */
308 arm_rfft_q15 (S->pRfft, pInlineBuffer, pState);
309
310 /*----------------------------------------------------------------------
311 * Step3: Multiply the FFT output with the weights.
312 *----------------------------------------------------------------------*/
313 arm_cmplx_mult_cmplx_q15 (pState, weights, pState, S->N);
314
315 /* The output of complex multiplication is in 3.13 format.
316 * Hence changing the format of N (i.e. 2*N elements) complex numbers to 1.15 format by shifting left by 2 bits. */
317 arm_shift_q15 (pState, 2, pState, S->N * 2);
318
319 /* ----------- Post-processing ---------- */
320 /* DCT-IV can be obtained from DCT-II by the equation,
321 * Y4(k) = Y2(k) - Y4(k-1) and Y4(-1) = Y4(0)
322 * Hence, Y4(0) = Y2(0)/2 */
323 /* Getting only real part from the output and Converting to DCT-IV */
324
325 /* pbuff initialized to input buffer. */
326 pbuff = pInlineBuffer;
327
328 /* pS1 initialized to pState */
329 pS1 = pState;
330
331 /* Calculating Y4(0) from Y2(0) using Y4(0) = Y2(0)/2 */
332 in = *pS1++ >> 1U;
333 /* input buffer acts as inplace, so output values are stored in the input itself. */
334 *pbuff++ = in;
335
336 /* pState pointer is incremented twice as the real values are located alternatively in the array */
337 pS1++;
338
339 /* Initializing the loop counter */
340 i = (S->N - 1U);
341
342 do
343 {
344 /* Calculating Y4(1) to Y4(N-1) from Y2 using equation Y4(k) = Y2(k) - Y4(k-1) */
345 /* pState pointer (pS1) is incremented twice as the real values are located alternatively in the array */
346 in = *pS1++ - in;
347 *pbuff++ = in;
348
349 /* points to the next real value */
350 pS1++;
351
352 /* Decrement loop counter */
353 i--;
354 } while (i > 0U);
355
356 /*------------ Normalizing the output by multiplying with the normalizing factor ----------*/
357
358 /* Initializing loop counter */
359 i = S->N;
360
361 /* pbuff initialized to the pInlineBuffer (now contains the output values) */
362 pbuff = pInlineBuffer;
363
364 do
365 {
366 /* Multiplying pInlineBuffer with the normalizing factor sqrt(2/N) */
367 in = *pbuff;
368 *pbuff++ = ((q15_t) (((q31_t) in * S->normalize) >> 15));
369
370 /* Decrement loop counter */
371 i--;
372
373 } while (i > 0U);
374
375 #endif /* #if defined (ARM_MATH_LOOPUNROLL) */
376
377 }
378
379 /**
380 @} end of DCT4_IDCT4 group
381 */
382