Mercurial > hg > ltpda

comparison m-toolbox/m/mdcs/mdc1_UTN/ltpda_series_derivative.m @ 0:f0afece42f48

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Import.
author Daniele Nicolodi <nicolodi@science.unitn.it>
date Wed, 23 Nov 2011 19:22:13 +0100
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:f0afece42f48
1 function varargout = ltpda_series_derivative(varargin)
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 %
4 % DESCRIPTION: LTPDA_SERIES_DERIVATIVE performs the numeric time derivative
5 % using the method of five points stencil Taylor series expansion [1].
6 % The function can perform the first, second, third and fourth derivetive
7 % of a series of input data.
8 %
9 %
10 %
11 % CALL: Deriv = ltpda_series_derivative(a, ... , pl)
12 %
13 %
14 % INPUTS:
15 % - Aos containing the data to be differentiated.
16 % - pl: is a parameter list containing the options for the
17 % calculation. Parameters are:
18 % - 'ORDER' whose accepeted options are:
19 % - 'FIRST' perform the first derivative using the
20 % couefficients vector d1 = [1 -8 0 8 -1]./(12*T)
21 % - 'SECOND' perform the second derivative using the
22 % coefficients vector d2 = [-1 16 -30 16 -1]./(12*T^2)
23 % - 'THIRD' perform the third derivative using the
24 % coefficients vector d3 = [-1 2 0 -2 1]./(2*T^3)
25 % - 'FOURTH' perform the third derivative using the
26 % coefficients vector d4 = [1 -4 6 -4 1]./(T^4)
27 %
28 % NOTE1: T is the sampling period
29 % NOTE2: The default option for 'ORDER' is 'SECOND'
30 %
31 %
32 %
33 % OUTPUTS:
34 % - Deriv is avector containing the resulting Aos after
35 % differentiation
36 %
37 %
38 % REFERENCES:
39 % [1] S. E. Koonin and D. C. Meredith, COMPUTATIONAL PHYSICS -
40 % Fortran Version, 1990, Westview Press
41 %
42 %
43 % VERSION: $Id: ltpda_series_derivative.m,v 1.1 2008/04/24 16:13:12 luigi Exp $
44 %
45 % HISTORY: 18-03-2008 L Ferraioli
46 % Creation
47 %
48 % The following call returns a parameter list object that contains the
49 % default parameter values:
50 %
51 % >> pl = ltpda_series_derivative('Params')
52 %
53 % The following call returns a string that contains the routine CVS version:
54 %
55 % >> version = ltpda_series_derivative('Version')
56 %
57 % The following call returns a string that contains the routine category:
58 %
59 % >> category = ltpda_series_derivative('Category')
60 %
61 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
62
63
64
65 %% Standard history variables
66
67 ALGONAME = mfilename;
68 VERSION = '$Id: ltpda_series_derivative.m,v 1.1 2008/04/24 16:13:12 luigi Exp $';
69 CATEGORY = 'Signal Processing';
70
71 %% Check if this is a call for parameters, the CVS version string
72 % or the method category
73 if nargin == 1 && ischar(varargin{1})
74 in = char(varargin{1});
75 if strcmpi(in, 'Params')
76 varargout{1} = getDefaultPL();
77 return
78 elseif strcmpi(in, 'Version')
79 varargout{1} = VERSION;
80 return
81 elseif strcmpi(in, 'Category')
82 varargout{1} = CATEGORY;
83 return
84 end
85 end
86
87 %% Collect input ao's, plist's and ao variable names
88
89 in_names = {};
90
91 for ii = 1:nargin
92 in_names{end+1} = inputname(ii); % Collect the inputs names
93 end
94
95 [as, pl, invars] = collect_inputs(varargin, in_names);
96
97 % produce one parameter list
98 if isa(pl, 'plist')
99 pl = combine(pl, getDefaultPL());
100 else
101 pl = getDefaultPL();
102 end
103
104 %% Initialize outputs
105 bs = [];
106
107 %% Go through analysis objects
108
109 % Assigning coefficients values based on the input options
110 switch find(pl, 'ORDER')
111 case 'FIRST'
112 Coeffs = [1 -8 0 8 -1]./12;
113 case 'SECOND'
114 Coeffs = [-1 16 -30 16 -1]./12;
115 case 'THIRD'
116 Coeffs = [-1 2 0 -2 1]./2;
117 case 'FOURTH'
118 Coeffs = [1 -4 6 -4 1];
119 end
120
121 for jj = 1:numel(as)
122
123 a = as(jj);
124 phi = a.data.y;
125 fs = get(a, 'fs');
126 T = 1/fs;
127
128 %% Building vectors for calculation
129 phi_temp = [phi(1);phi(1);phi(1);phi(1);phi;phi(end);phi(end);phi(end);phi(end)];
130
131 % Switching between the input options
132 switch find(pl, 'ORDER')
133 case 'FIRST'
134 Deriv = (1/T).*(Coeffs(1)*phi_temp(1:end-4) + Coeffs(2)*phi_temp(2:end-3) + Coeffs(3)*phi_temp(3:end-2) + Coeffs(4)*phi_temp(4:end-1) + Coeffs(5)*phi_temp(5:end));
135 Deriv = Deriv(3:end-2);
136 case 'SECOND'
137 Deriv = (1/T^2).*(Coeffs(1)*phi_temp(1:end-4) + Coeffs(2)*phi_temp(2:end-3) + Coeffs(3)*phi_temp(3:end-2) + Coeffs(4)*phi_temp(4:end-1) + Coeffs(5)*phi_temp(5:end));
138 Deriv = Deriv(3:end-2);
139 case 'THIRD'
140 Deriv = (1/T^3).*(Coeffs(1)*phi_temp(1:end-4) + Coeffs(2)*phi_temp(2:end-3) + Coeffs(3)*phi_temp(3:end-2) + Coeffs(4)*phi_temp(4:end-1) + Coeffs(5)*phi_temp(5:end));
141 Deriv = Deriv(3:end-2);
142 case 'FOURTH'
143 Deriv = (1/T^4).*(Coeffs(1)*phi_temp(1:end-4) + Coeffs(2)*phi_temp(2:end-3) + Coeffs(3)*phi_temp(3:end-2) + Coeffs(4)*phi_temp(4:end-1) + Coeffs(5)*phi_temp(5:end));
144 Deriv = Deriv(3:end-2);
145 end
146
147 %% Building output Aos
148
149 % Creating new time series data
150 Deriv = tsdata(Deriv,fs);
151
152 Deriv = set(Deriv, 'name', ...
153 sprintf('ltpda_parfit_derivative(%s)', a.data.name), ...
154 'xunits', a.data.xunits, 'yunits', 'm/s^2');
155
156 % make new history objects
157 h = history(ALGONAME, VERSION, pl, a.hist);
158 h = set(h, 'invars', invars);
159
160 % make output analysis objects
161 b = ao(Deriv, h);
162 clear Deriv
163
164 % name, mfilename, description for these objects
165 b = setnh(b, 'name', sprintf('ltpda_parfit_derivative(%s)', invars{jj}),...
166 'description', find(pl,'description'));
167
168 %% Output the data
169
170 bs = [bs b];
171 clear a b
172
173 end
174
175 %% Reshape the ouput to the same size of the input
176 varargout{1} = reshape(bs, size(as));
177
178 %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
179 %
180 % FUNCTION: getDefaultPL
181 %
182 % DESCRIPTION: Get default params
183 %
184 % HISTORY: 01-02-2008 M Hueller
185 % Creation
186 %
187 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
188 function plo = getDefaultPL()
189
190 plo = plist('NAME', 'ltpda_parfit_derivative', 'ORDER', 'SECOND', ...
191 'DESCRIPTION','Series Expansion Derivative Method');
192
193
194 % END
195
196
197
198