libdspl-2.0
Digital Signal Processing Algorithm Library
statistic.c
1 /*
2 * Copyright (c) 2015-2019 Sergey Bakhurin
3 * Digital Signal Processing Library [http://dsplib.org]
4 *
5 * This file is part of DSPL.
6 *
7 * is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * DSPL is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with Foobar. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #include <stdio.h>
23 #include <stdlib.h>
24 #include <string.h>
25 #include "dspl.h"
26 
27 
28 
29 
30 
31 #ifdef DOXYGEN_ENGLISH
32 
77 #endif
78 #ifdef DOXYGEN_RUSSIAN
79 
122 #endif
123 int DSPL_API find_max_abs(double* a, int n, double* m, int* ind)
124 {
125  int k, i;
126  double t;
127  if(!a)
128  return ERROR_PTR;
129  if(n < 1)
130  return ERROR_SIZE;
131  t = fabs(a[0]);
132  i = 0;
133  for(k = 1; k < n; k++)
134  {
135  if(fabs(a[k]) > t)
136  {
137  t = fabs(a[k]);
138  i = k;
139  }
140  }
141  if(m)
142  *m = t;
143  if(ind)
144  *ind = i;
145  return RES_OK;
146 }
147 
148 
149 #ifdef DOXYGEN_ENGLISH
150 
151 #endif
152 #ifdef DOXYGEN_RUSSIAN
153 
154 #endif
155 int DSPL_API histogram(double* x, int n, int nh, double* pedges, double* ph)
156 {
157  double xmin, xmax;
158  int k, ind;
159  int res;
160 
161  if(!x || !pedges || !ph)
162  return ERROR_PTR;
163 
164  if(n<1 || nh<1)
165  return ERROR_SIZE;
166 
167  res = minmax(x, n, &xmin, &xmax);
168  if(res != RES_OK)
169  return res;
170 
171  res = linspace(xmin, xmax, nh+1, DSPL_SYMMETRIC, pedges);
172  if(res != RES_OK)
173  return res;
174 
175  memset(ph, 0, nh*sizeof(double));
176  for(k = 0; k < n; k++)
177  {
178  ind = 0;
179  while(ind<nh && x[k]>=pedges[ind])
180  ind++;
181  ph[ind-1]+=1.0;
182  }
183  return RES_OK;
184 }
185 
186 
187 
188 
189 #ifdef DOXYGEN_ENGLISH
190 
191 #endif
192 #ifdef DOXYGEN_RUSSIAN
193 
194 #endif
195 int DSPL_API histogram_norm(double* y, int n, int nh, double* x, double* w)
196 {
197  double *pedges = NULL;
198  int k, res;
199 
200  if(!y || !x || !w)
201  return ERROR_PTR;
202 
203  if(n<1 || nh<1)
204  return ERROR_SIZE;
205 
206  pedges = (double*)malloc((nh+1)*sizeof(double));
207 
208  res = histogram(y, n, nh, pedges, w);
209  if(res != RES_OK)
210  goto exit_label;
211 
212  for(k = 1; k < nh+1; k++)
213  {
214  x[k-1] = 0.5*(pedges[k] + pedges[k-1]);
215  w[k-1] /= ((double)n * (pedges[k] - pedges[k-1]));
216  }
217 
218  res = RES_OK;
219 exit_label:
220 
221  if(pedges)
222  free(pedges);
223  return res;
224 
225 }
226 
227 
228 
229 
230 
231 
232 
233 #ifdef DOXYGEN_ENGLISH
234 
235 #endif
236 #ifdef DOXYGEN_RUSSIAN
237 
238 #endif
239 int DSPL_API mean(double* x, int n, double* m)
240 {
241  int k;
242  double sum;
243  if(!x || !m)
244  return ERROR_PTR;
245  if(n<1)
246  return ERROR_SIZE;
247 
248  sum = x[0];
249  for(k = 1; k < n; k++)
250  sum += x[k];
251 
252  (*m) = sum / (double)n;
253  return RES_OK;
254 }
255 
256 
257 
258 
259 
260 #ifdef DOXYGEN_ENGLISH
261 
262 #endif
263 #ifdef DOXYGEN_RUSSIAN
264 
265 #endif
266 int DSPL_API mean_cmplx(complex_t* x, int n, complex_t* m)
267 {
268  int k;
269  complex_t sum;
270  if(!x || !m)
271  return ERROR_PTR;
272  if(n<1)
273  return ERROR_SIZE;
274 
275  RE(sum) = RE(x[0]);
276  IM(sum) = IM(x[0]);
277  for(k = 1; k < n; k++)
278  {
279  RE(sum) += RE(x[k]);
280  IM(sum) += IM(x[k]);
281  }
282  RE(sum) /= (double)n;
283  IM(sum) /= (double)n;
284  memcpy(m, sum, sizeof(complex_t));
285 
286  return RES_OK;
287 }
288 
289 
290 
291 
292 #ifdef DOXYGEN_ENGLISH
293 
294 #endif
295 #ifdef DOXYGEN_RUSSIAN
296 
297 #endif
298 int DSPL_API minmax(double* x, int n, double* xmin, double* xmax)
299 {
300  int k;
301  double min, max;
302 
303  if(!x)
304  return ERROR_PTR;
305 
306  if(n<1)
307  return ERROR_SIZE;
308 
309  min = max = x[0];
310  for(k = 1; k < n; k++)
311  {
312  min = x[k] < min ? x[k] : min;
313  max = x[k] > max ? x[k] : max;
314  }
315 
316  if(xmin)
317  *xmin = min;
318  if(xmax)
319  *xmax = max;
320 
321  return RES_OK;
322 }
323 
324 
325 #ifdef DOXYGEN_ENGLISH
326 
327 #endif
328 #ifdef DOXYGEN_RUSSIAN
329 
330 #endif
331 int DSPL_API std(double* x, int n, double* s)
332 {
333  int k, err;
334  double sum, m;
335  err = mean(x, n, &m);
336  if(err != RES_OK)
337  goto exit_label;
338  sum = (x[0] - m) * (x[0] - m);
339  for(k = 1; k < n; k++)
340  sum += (x[k] - m) * (x[k] - m);
341  (*s) = sqrt(sum / (double)(n-1));
342 exit_label:
343  return err;
344 }
345 
346 
347 
348 
349 
350 #ifdef DOXYGEN_ENGLISH
351 
352 #endif
353 #ifdef DOXYGEN_RUSSIAN
354 
355 #endif
356 int DSPL_API std_cmplx(complex_t* x, int n, double* s)
357 {
358  int k, err;
359  complex_t tmp, m;
360  double sum;
361  err = mean_cmplx(x, n, &m);
362  if(err != RES_OK)
363  goto exit_label;
364 
365  RE(tmp) = RE(x[0]) - RE(m);
366  IM(tmp) = IM(x[0]) - IM(m);
367  sum = ABSSQR(tmp);
368  for(k = 1; k < n; k++)
369  {
370  RE(tmp) = RE(x[k]) - RE(m);
371  IM(tmp) = IM(x[k]) - IM(m);
372  sum += ABSSQR(tmp);
373  }
374  *s = sqrt(sum / (double)(n-1));
375 exit_label:
376  return err;
377 }
378 
379 
380 
381 
#define RE(x)
Macro sets real part of the complex number.
Definition: dspl.h:359
#define ERROR_PTR
Pointer error. This error means that one of the required pointers (memory to be allocated for) is tra...
Definition: dspl.h:545
#define ERROR_SIZE
Error array size. This error occurs when in addition to the pointer the wrong input is passed to the ...
Definition: dspl.h:553
double complex_t[2]
Complex data type.
Definition: dspl.h:86
int DSPL_API find_max_abs(double *a, int n, double *m, int *ind)
Find maximum absolute value from the real vector a
Definition: statistic.c:123
#define RES_OK
The function completed correctly. No errors.
Definition: dspl.h:497
#define ABSSQR(x)
The macro returns the square of the modulus of a complex number x.
Definition: dspl.h:473
int DSPL_API linspace(double x0, double x1, int n, int type, double *x)
Function fills a vector with n linearly spaced elements between x0 and x1.
Definition: array.c:1009
#define IM(x)
Macro sets imaginary part of the complex number.
Definition: dspl.h:417