Merge branch 'jnosky/master'
[fw/stlink] / exampleF4 / CMSIS / DSP_Lib / Source / StatisticsFunctions / arm_var_q31.c
1 /* ----------------------------------------------------------------------   
2 * Copyright (C) 2010 ARM Limited. All rights reserved.   
3 *   
4 * $Date:        15. July 2011  
5 * $Revision:    V1.0.10  
6 *   
7 * Project:          CMSIS DSP Library   
8 * Title:                arm_var_q31.c   
9 *   
10 * Description:  Variance of an array of Q31 type.   
11 *   
12 * Target Processor: Cortex-M4/Cortex-M3/Cortex-M0
13 *  
14 * Version 1.0.10 2011/7/15 
15 *    Big Endian support added and Merged M0 and M3/M4 Source code.  
16 *   
17 * Version 1.0.3 2010/11/29  
18 *    Re-organized the CMSIS folders and updated documentation.   
19 *    
20 * Version 1.0.2 2010/11/11   
21 *    Documentation updated.    
22 *   
23 * Version 1.0.1 2010/10/05    
24 *    Production release and review comments incorporated.   
25 *   
26 * Version 1.0.0 2010/09/20    
27 *    Production release and review comments incorporated.   
28 * -------------------------------------------------------------------- */
29
30 #include "arm_math.h"
31
32 /**   
33  * @ingroup groupStats   
34  */
35
36 /**   
37  * @addtogroup variance   
38  * @{   
39  */
40
41 /**   
42  * @brief Variance of the elements of a Q31 vector.   
43  * @param[in]       *pSrc points to the input vector   
44  * @param[in]       blockSize length of the input vector   
45  * @param[out]      *pResult variance value returned here   
46  * @return none.   
47  *   
48  * @details   
49  * <b>Scaling and Overflow Behavior:</b>   
50  *   
51  *\par   
52  * The function is implemented using an internal 64-bit accumulator.   
53  * The input is represented in 1.31 format, and intermediate multiplication   
54  * yields a 2.62 format.   
55  * The accumulator maintains full precision of the intermediate multiplication results,    
56  * but provides only a single guard bit.   
57  * There is no saturation on intermediate additions.   
58  * If the accumulator overflows it wraps around and distorts the result.   
59  * In order to avoid overflows completely the input signal must be scaled down by    
60  * log2(blockSize) bits, as a total of blockSize additions are performed internally.    
61  * Finally, the 2.62 accumulator is right shifted by 31 bits to yield a 1.31 format value.   
62  *   
63  */
64
65
66 void arm_var_q31(
67   q31_t * pSrc,
68   uint32_t blockSize,
69   q63_t * pResult)
70 {
71   q63_t sum = 0;                                 /* Accumulator */
72   q31_t meanOfSquares, squareOfMean;             /* Mean of square and square of mean */
73   q31_t mean;                                    /* Mean */
74   q31_t in;                                      /* Input variable */
75   q31_t t;                                       /* Temporary variable */
76   uint32_t blkCnt;                               /* loop counter */
77
78 #ifndef ARM_MATH_CM0
79
80   /* Run the below code for Cortex-M4 and Cortex-M3 */
81
82   q31_t *pIn;                                    /* Temporary pointer */
83
84   pIn = pSrc;
85
86   /*loop Unrolling */
87   blkCnt = blockSize >> 2u;
88
89   /* First part of the processing with loop unrolling.  Compute 4 outputs at a time.   
90    ** a second loop below computes the remaining 1 to 3 samples. */
91   while(blkCnt > 0u)
92   {
93     /* C = (A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1])  */
94     /* Compute Sum of squares of the input samples   
95      * and then store the result in a temporary variable, sum. */
96     in = *pSrc++;
97     sum += ((q63_t) (in) * (in));
98     in = *pSrc++;
99     sum += ((q63_t) (in) * (in));
100     in = *pSrc++;
101     sum += ((q63_t) (in) * (in));
102     in = *pSrc++;
103     sum += ((q63_t) (in) * (in));
104
105     /* Decrement the loop counter */
106     blkCnt--;
107   }
108
109   /* If the blockSize is not a multiple of 4, compute any remaining output samples here.   
110    ** No loop unrolling is used. */
111   blkCnt = blockSize % 0x4u;
112
113   while(blkCnt > 0u)
114   {
115     /* C = (A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1]) */
116     /* Compute Sum of squares of the input samples   
117      * and then store the result in a temporary variable, sum. */
118     in = *pSrc++;
119     sum += ((q63_t) (in) * (in));
120
121     /* Decrement the loop counter */
122     blkCnt--;
123   }
124
125   /* Compute Mean of squares of the input samples   
126    * and then store the result in a temporary variable, meanOfSquares. */
127   t = (q31_t) ((1.0 / (blockSize - 1)) * 1073741824LL);
128   sum = (sum >> 31);
129   meanOfSquares = (q31_t) ((sum * t) >> 30);
130
131   /* Reset the accumulator */
132   sum = 0;
133
134   /*loop Unrolling */
135   blkCnt = blockSize >> 2u;
136
137   /* Reset the input working pointer */
138   pSrc = pIn;
139
140   /* First part of the processing with loop unrolling.  Compute 4 outputs at a time.   
141    ** a second loop below computes the remaining 1 to 3 samples. */
142   while(blkCnt > 0u)
143   {
144     /* C = (A[0] + A[1] + A[2] + ... + A[blockSize-1]) */
145     /* Compute sum of all input values and then store the result in a temporary variable, sum. */
146     sum += *pSrc++;
147     sum += *pSrc++;
148     sum += *pSrc++;
149     sum += *pSrc++;
150
151     /* Decrement the loop counter */
152     blkCnt--;
153   }
154
155   /* If the blockSize is not a multiple of 4, compute any remaining output samples here.   
156    ** No loop unrolling is used. */
157   blkCnt = blockSize % 0x4u;
158
159   while(blkCnt > 0u)
160   {
161     /* C = (A[0] + A[1] + A[2] + ... + A[blockSize-1]) */
162     /* Compute sum of all input values and then store the result in a temporary variable, sum. */
163     sum += *pSrc++;
164
165     /* Decrement the loop counter */
166     blkCnt--;
167   }
168
169 #else
170
171   /* Run the below code for Cortex-M0 */
172
173   q63_t sumOfSquares = 0;                        /* Accumulator */
174   /* Loop over blockSize number of values */
175   blkCnt = blockSize;
176
177   while(blkCnt > 0u)
178   {
179     /* C = (A[0] * A[0] + A[1] * A[1] + ... + A[blockSize-1] * A[blockSize-1])  */
180     /* Compute Sum of squares of the input samples    
181      * and then store the result in a temporary variable, sumOfSquares. */
182     in = *pSrc++;
183     sumOfSquares += ((q63_t) (in) * (in));
184
185     /* C = (A[0] + A[1] + A[2] + ... + A[blockSize-1]) */
186     /* Compute sum of all input values and then store the result in a temporary variable, sum. */
187     sum += in;
188
189     /* Decrement the loop counter */
190     blkCnt--;
191   }
192
193   /* Compute Mean of squares of the input samples    
194    * and then store the result in a temporary variable, meanOfSquares. */
195   t = (q31_t) ((1.0 / (blockSize - 1)) * 1073741824LL);
196   sumOfSquares = (sumOfSquares >> 31);
197   meanOfSquares = (q31_t) ((sumOfSquares * t) >> 30);
198
199 #endif /* #ifndef ARM_MATH_CM0 */
200
201   /* Compute mean of all input values */
202   t = (q31_t) ((1.0 / (blockSize * (blockSize - 1u))) * 2147483648LL);
203   mean = (q31_t) (sum);
204
205   /* Compute square of mean */
206   squareOfMean = (q31_t) (((q63_t) mean * mean) >> 31);
207   squareOfMean = (q31_t) (((q63_t) squareOfMean * t) >> 31);
208
209   /* Compute variance and then store the result to the destination */
210   *pResult = (q63_t) meanOfSquares - squareOfMean;
211
212 }
213
214 /**   
215  * @} end of variance group   
216  */