1 /* ----------------------------------------------------------------------
2 * Project: CMSIS DSP Library
3 * Title: arm_std_q15.c
4 * Description: Standard deviation of an array of Q15 vector
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/statistics_functions.h"
30
31 /**
32 @ingroup groupStats
33 */
34
35 /**
36 @addtogroup STD
37 @{
38 */
39
40 /**
41 @brief Standard deviation of the elements of a Q15 vector.
42 @param[in] pSrc points to the input vector
43 @param[in] blockSize number of samples in input vector
44 @param[out] pResult standard deviation value returned here
45 @return none
46
47 @par Scaling and Overflow Behavior
48 The function is implemented using a 64-bit internal accumulator.
49 The input is represented in 1.15 format.
50 Intermediate multiplication yields a 2.30 format, and this
51 result is added without saturation to a 64-bit accumulator in 34.30 format.
52 With 33 guard bits in the accumulator, there is no risk of overflow, and the
53 full precision of the intermediate multiplication is preserved.
54 Finally, the 34.30 result is truncated to 34.15 format by discarding the lower
55 15 bits, and then saturated to yield a result in 1.15 format.
56 */
57 #if defined(ARM_MATH_MVEI) && !defined(ARM_MATH_AUTOVECTORIZE)
arm_std_q15(const q15_t * pSrc,uint32_t blockSize,q15_t * pResult)58 void arm_std_q15(
59 const q15_t * pSrc,
60 uint32_t blockSize,
61 q15_t * pResult)
62 {
63 q15_t var=0;
64
65 arm_var_q15(pSrc, blockSize, &var);
66 arm_sqrt_q15(var,pResult);
67 }
68 #else
arm_std_q15(const q15_t * pSrc,uint32_t blockSize,q15_t * pResult)69 void arm_std_q15(
70 const q15_t * pSrc,
71 uint32_t blockSize,
72 q15_t * pResult)
73 {
74 uint32_t blkCnt; /* Loop counter */
75 q31_t sum = 0; /* Accumulator */
76 q31_t meanOfSquares, squareOfMean; /* Square of mean and mean of square */
77 q63_t sumOfSquares = 0; /* Sum of squares */
78 q15_t in; /* Temporary variable to store input value */
79
80 #if defined (ARM_MATH_LOOPUNROLL) && defined (ARM_MATH_DSP)
81 q31_t in32; /* Temporary variable to store input value */
82 #endif
83
84 if (blockSize <= 1U)
85 {
86 *pResult = 0;
87 return;
88 }
89
90 #if defined (ARM_MATH_LOOPUNROLL)
91
92 /* Loop unrolling: Compute 4 outputs at a time */
93 blkCnt = blockSize >> 2U;
94
95 while (blkCnt > 0U)
96 {
97 /* C = A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1] */
98 /* C = A[0] + A[1] + ... + A[blockSize-1] */
99
100 /* Compute sum of squares and store result in a temporary variable, sumOfSquares. */
101 /* Compute sum and store result in a temporary variable, sum. */
102 #if defined (ARM_MATH_DSP)
103 in32 = read_q15x2_ia ((q15_t **) &pSrc);
104 sumOfSquares = __SMLALD(in32, in32, sumOfSquares);
105 sum += ((in32 << 16U) >> 16U);
106 sum += (in32 >> 16U);
107
108 in32 = read_q15x2_ia ((q15_t **) &pSrc);
109 sumOfSquares = __SMLALD(in32, in32, sumOfSquares);
110 sum += ((in32 << 16U) >> 16U);
111 sum += (in32 >> 16U);
112 #else
113 in = *pSrc++;
114 sumOfSquares += (in * in);
115 sum += in;
116
117 in = *pSrc++;
118 sumOfSquares += (in * in);
119 sum += in;
120
121 in = *pSrc++;
122 sumOfSquares += (in * in);
123 sum += in;
124
125 in = *pSrc++;
126 sumOfSquares += (in * in);
127 sum += in;
128 #endif /* #if defined (ARM_MATH_DSP) */
129
130 /* Decrement loop counter */
131 blkCnt--;
132 }
133
134 /* Loop unrolling: Compute remaining outputs */
135 blkCnt = blockSize % 0x4U;
136
137 #else
138
139 /* Initialize blkCnt with number of samples */
140 blkCnt = blockSize;
141
142 #endif /* #if defined (ARM_MATH_LOOPUNROLL) */
143
144 while (blkCnt > 0U)
145 {
146 /* C = A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1] */
147 /* C = A[0] + A[1] + ... + A[blockSize-1] */
148
149 in = *pSrc++;
150 /* Compute sum of squares and store result in a temporary variable, sumOfSquares. */
151 sumOfSquares += (in * in);
152 /* Compute sum and store result in a temporary variable, sum. */
153 sum += in;
154
155 /* Decrement loop counter */
156 blkCnt--;
157 }
158
159 /* Compute Mean of squares and store result in a temporary variable, meanOfSquares. */
160 meanOfSquares = (q31_t) (sumOfSquares / (q63_t)(blockSize - 1U));
161
162 /* Compute square of mean */
163 squareOfMean = (q31_t) ((q63_t) sum * sum / (q63_t)(blockSize * (blockSize - 1U)));
164
165 /* mean of squares minus the square of mean. */
166 /* Compute standard deviation and store result in destination */
167 arm_sqrt_q15(__SSAT((meanOfSquares - squareOfMean) >> 15U, 16U), pResult);
168 }
169 #endif /* defined(ARM_MATH_MVEI) */
170
171 /**
172 @} end of STD group
173 */
174
