File indexing completed on 2024-05-12 03:47:54
0001 /* 0002 File : nsl_smooth.h 0003 Project : LabPlot 0004 Description : NSL smooth functions 0005 -------------------------------------------------------------------- 0006 SPDX-FileCopyrightText: 2016 Stefan Gerlach <stefan.gerlach@uni.kn> 0007 SPDX-License-Identifier: GPL-2.0-or-later 0008 */ 0009 0010 #ifndef NSL_SMOOTH_H 0011 #define NSL_SMOOTH_H 0012 0013 #include <gsl/gsl_matrix.h> 0014 0015 #define NSL_SMOOTH_TYPE_COUNT 4 0016 typedef enum { 0017 nsl_smooth_type_moving_average, 0018 nsl_smooth_type_moving_average_lagged, 0019 nsl_smooth_type_percentile, 0020 nsl_smooth_type_savitzky_golay 0021 } nsl_smooth_type; 0022 extern const char* nsl_smooth_type_name[]; 0023 /* TODO: LOWESS/etc., Bezier, B-Spline, (FFT Filter) */ 0024 0025 /* mode of extension for padding signal 0026 * none: reduce points at edges 0027 * interp: polynomial interpolation 0028 * mirror: 3 2 | 1 2 3 4 5 | 4 3 (reflect) 0029 * nearest: 1 1 | 1 2 3 4 5 | 5 5 (repeat) 0030 * constant: L L | 1 2 3 4 5 | R R 0031 * periodic: 4 5 | 1 2 3 4 5 | 1 2 (wrap) 0032 */ 0033 #define NSL_SMOOTH_PAD_MODE_COUNT 6 0034 typedef enum { 0035 nsl_smooth_pad_none, 0036 nsl_smooth_pad_interp, 0037 nsl_smooth_pad_mirror, 0038 nsl_smooth_pad_nearest, 0039 nsl_smooth_pad_constant, 0040 nsl_smooth_pad_periodic 0041 } nsl_smooth_pad_mode; 0042 extern const char* nsl_smooth_pad_mode_name[]; 0043 0044 #define NSL_SMOOTH_WEIGHT_TYPE_COUNT 8 0045 typedef enum { 0046 nsl_smooth_weight_uniform, 0047 nsl_smooth_weight_triangular, 0048 nsl_smooth_weight_binomial, 0049 nsl_smooth_weight_parabolic, 0050 nsl_smooth_weight_quartic, 0051 nsl_smooth_weight_triweight, 0052 nsl_smooth_weight_tricube, 0053 nsl_smooth_weight_cosine 0054 } nsl_smooth_weight_type; 0055 extern const char* nsl_smooth_weight_type_name[]; 0056 /*TODO: IIR: exponential, Gaussian, see nsl_sf_kernel */ 0057 0058 /* values used for constant padding */ 0059 extern double nsl_smooth_pad_constant_lvalue, nsl_smooth_pad_constant_rvalue; 0060 0061 /********* Smoothing algorithms **********/ 0062 0063 /* Moving average */ 0064 int nsl_smooth_moving_average(double* data, size_t n, size_t points, nsl_smooth_weight_type weight, nsl_smooth_pad_mode mode); 0065 0066 /* Lagged moving average */ 0067 int nsl_smooth_moving_average_lagged(double* data, size_t n, size_t points, nsl_smooth_weight_type weight, nsl_smooth_pad_mode mode); 0068 0069 /* Percentile filter */ 0070 int nsl_smooth_percentile(double* data, size_t n, size_t points, double percentile, nsl_smooth_pad_mode mode); 0071 0072 /* Savitzky-Golay coefficients */ 0073 /** 0074 * \brief Compute Savitzky-Golay coefficients and store them into #h. 0075 * 0076 * This function follows GSL conventions in that it writes its result into a matrix allocated by 0077 * the caller and returns a non-zero result on error. 0078 * 0079 * The coefficient matrix is defined as the matrix H mapping a set of input values to the values 0080 * of the polynomial of order #polynom_order which minimizes squared deviations from the input 0081 * values. It is computed using the formula \$H=V(V^TV)^(-1)V^T\$, where \$V\$ is the Vandermonde 0082 * matrix of the point indices. 0083 * 0084 * For a short description of the mathematical background, see 0085 * http://www.statistics4u.info/fundstat_eng/cc_filter_savgol_math.html 0086 */ 0087 int nsl_smooth_savgol_coeff(size_t points, int order, gsl_matrix* h); 0088 0089 /* set values for constant padding */ 0090 void nsl_smooth_pad_constant_set(double lvalue, double rvalue); 0091 0092 /* Savitzky-Golay smoothing */ 0093 /** 0094 * \brief Savitzky-Golay smoothing of (uniformly distributed) data. 0095 * 0096 * When the data is not uniformly distributed, Savitzky-Golay looses its interesting conservation 0097 * properties. On the other hand, a central point of the algorithm is that for uniform data, the 0098 * operation can be implemented as a convolution. This is considerably more efficient than a more 0099 * generic method able to handle non-uniform input data. 0100 */ 0101 int nsl_smooth_savgol(double* data, size_t n, size_t points, int order, nsl_smooth_pad_mode mode); 0102 0103 /* Savitzky-Golay default smoothing (interp) */ 0104 int nsl_smooth_savgol_default(double* data, size_t n, size_t points, int order); 0105 0106 /* TODO SmoothFilter::smoothModifiedSavGol(double *x_in, double *y_inout) 0107 see SmoothFilter.cpp of libscidavis 0108 */ 0109 0110 #endif /* NSL_SMOOTH_H */