1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
|
function M = cummax(varargin)
// Cumulative maximum
//
// Calling Sequence
// M = cummax(A)
// returns the cumulative maximum of the arguments of A. The dimension
// of M is same as the dimension of A. If A is a 2D matrix, the operation
// is performed along the columns. For a hypermatrix, the operation is
// performed along the first non-zero dimension
// M = cummax(A,dim)
// The operation is performed along the dimension specified by dim
// M = cummax(_,direction)
// direction specifies as the direction of operation
//
// Parameters
// A - real|complex numbers - vector|matrix
// Input Array
// For complex elements, cummax compares the magnitude of elements. If
// the magnitude are same, phase angles are compared.
// dim - positive integer - scalar
// Dimension to operate along
// If no dimension is specified, then the default value is the first
// array dimension whose value is greater than 1
// direction - string flag - 'forward' (default) or 'reverse'
// Direction of cumulation
// If the direction is forward, cummax works from 1 to end of the active
// dimension. Otherwise, it works in the opposite sense
//
// Examples
// 1) Cumulative maximum values in a vector
// v = [8 9 1 10 6 1 3 6 10 10]
// M = cummax(v)
//
// Expected output: [8 9 9 10 10 10 10 10 10 10]
//
// Authors
// Ayush Baid
//
// See Also
// cummax | cumprod | cumsum | max | max
[numOutArgs,numInArgs] = argn(0);
// ** Checking number of arguments
if numInArgs<1 | numInArgs>3 then
msg = "cummax: Wrong number of input argument; 1-6 expected";
error(77,msg);
end
if numOutArgs~=1 then
msg = "cummax: Wrong number of output argument; 1 expected";
error(78,msg);
end
// ** Parsing input args **
// defining default arguments
isForward = %t;
dim = [];
directionArg = "";
A = varargin(1);
// A should contain numeric entries
if ~(type(A)==1 | type(A)==8 | type(A)==17) then
msg = "cummax: Wrong type for argument #1 (A); Real or complex entries expected ";
error(53,msg);
end
if numInArgs>1 then
temp = varargin(2);
if type(temp)==10 then
// it is the direction argument
directionArg = temp;
elseif type(temp)==1 | type(temp)==8 then
dim = int(temp);
else
msg = "cummax: Wrong type for argument #2; Either dim (integer) or direction (string) expected";
error(53,msg);
end
end
if numInArgs>2 then
directionArg = varargin(3);
if type(directionArg)~=10 then
msg = "cummax: Wrong type for argument #3 (direction); String expected";
error(53,msg);
end
end
if isempty(dim) then
dimArray = 1:ndims(A);
dim = find(size(A)~=1,1);
end
// additional checks on dim
if size(A,dim)==1 then
M = A;
return
end
// extracting direction
if strcmpi(directionArg,"reverse")==0 then
isForward = %f;
elseif strcmpi(directionArg,"forward")==0 then
isForward = %t;
elseif strcmpi(directionArg,"")~=0 then
msg = "cummax: Wrong value for argument #3 (direction)";
error(53,msg);
end
sizeA = size(A);
sizeDim = size(A,dim);
// restructuring A into a 3D matrix with the specified dimension as the middle elements
leftSize = prod(sizeA(1:dim-1));
rightSize = prod(sizeA(dim+1:$));
middleSize = sizeDim;
A_ = matrix(A,[leftSize,middleSize,rightSize]);
M_ = zeros(leftSize,middleSize,rightSize);
for i=1:leftSize
for j=1:rightSize
M_(i,:,j) = cummaxVec(A_(i,:,j),isForward);
end
end
M = matrix(M_,sizeA);
endfunction
function out = cummaxVec(inp,isForward)
// performs cummax on vector inputs
if isForward then
startIndex=1;
endIndex = length(inp);
step = 1;
else
startIndex=length(inp);
endIndex = 1;
step = -1;
end
out(startIndex) = inp(startIndex);
if isreal(inp) then
for i=startIndex+step:step:endIndex
if isnan(out(i-step)) then
out(i) = inp(i);
elseif inp(i)>=out(i-step) then
out(i) = inp(i);
else
out(i) = out(i-step);
end
end
else
magVec = abs(inp);
phaseVec = atan(imag(inp),real(inp));
// phase - first compare absolute value; then give priority to positive phases
prevMag = magVec(startIndex);
prevPhase = phaseVec(startIndex);
for i=(startIndex+step):step:endIndex
if isnan(out(i-step)) then
out(i) = inp(i);
prevMag = magVec(i);
prevPhase = phaseVec(i);
elseif magVec(i)>prevMag then
out(i) = inp(i);
prevMag = magVec(i);
prevPhase = phaseVec(i);
elseif magVec(i)<prevMag then
out(i) = out(i-step);
else
if phaseVec(i)>prevPhase then
out(i) = inp(i);
prevMag = magVec(i);
prevPhase = phaseVec(i);
else
out(i) = out(i-step);
end
end
end
end
endfunction
|
