|
2928
|
1 /* |
|
|
2 |
|
|
3 Copyright (C) 1996, 1997 John W. Eaton |
|
4850
|
4 Copyright (C) 2004 David Bateman |
|
2928
|
5 |
|
|
6 This file is part of Octave. |
|
|
7 |
|
|
8 Octave is free software; you can redistribute it and/or modify it |
|
|
9 under the terms of the GNU General Public License as published by the |
|
|
10 Free Software Foundation; either version 2, or (at your option) any |
|
|
11 later version. |
|
|
12 |
|
|
13 Octave is distributed in the hope that it will be useful, but WITHOUT |
|
|
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
|
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
|
16 for more details. |
|
|
17 |
|
|
18 You should have received a copy of the GNU General Public License |
|
|
19 along with Octave; see the file COPYING. If not, write to the Free |
|
|
20 Software Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. |
|
|
21 |
|
|
22 */ |
|
|
23 |
|
|
24 #ifdef HAVE_CONFIG_H |
|
|
25 #include <config.h> |
|
|
26 #endif |
|
|
27 |
|
3826
|
28 #include "lo-mappers.h" |
|
4153
|
29 #include "quit.h" |
|
3826
|
30 |
|
2928
|
31 #include "defun-dld.h" |
|
|
32 #include "error.h" |
|
|
33 #include "gripes.h" |
|
|
34 #include "oct-obj.h" |
|
4850
|
35 #include "lo-ieee.h" |
|
|
36 #include "data-conv.h" |
|
|
37 #include "ov-cx-mat.h" |
|
|
38 #include "oct-sort.cc" |
|
2928
|
39 |
|
4853
|
40 // If we have IEEE 754 data format, then we can use the trick of |
|
|
41 // casting doubles as unsigned eight byte integers, and with a little |
|
|
42 // bit of magic we can automatically sort the NaN's correctly. |
|
2928
|
43 |
|
4853
|
44 #if defined (HAVE_IEEE754_DATA_FORMAT) && defined (EIGHT_BYTE_INT) |
|
4850
|
45 |
|
|
46 static inline unsigned EIGHT_BYTE_INT FloatFlip(unsigned EIGHT_BYTE_INT f) |
|
2928
|
47 { |
|
4850
|
48 unsigned EIGHT_BYTE_INT mask = -(EIGHT_BYTE_INT)(f >> 63) | |
|
|
49 0x8000000000000000ULL; |
|
|
50 return f ^ mask; |
|
|
51 } |
|
2928
|
52 |
|
4850
|
53 inline unsigned EIGHT_BYTE_INT IFloatFlip(unsigned EIGHT_BYTE_INT f) |
|
|
54 { |
|
|
55 unsigned EIGHT_BYTE_INT mask = ((f >> 63) - 1) | 0x8000000000000000ULL; |
|
|
56 return f ^ mask; |
|
|
57 } |
|
2928
|
58 |
|
4850
|
59 struct vec_index |
|
|
60 { |
|
|
61 unsigned EIGHT_BYTE_INT vec; |
|
|
62 int indx; |
|
|
63 }; |
|
2928
|
64 |
|
4850
|
65 bool |
|
|
66 ieee754_compare (vec_index *a, vec_index *b) |
|
|
67 { |
|
|
68 return (a->vec < b->vec); |
|
2928
|
69 } |
|
|
70 |
|
4883
|
71 template class octave_sort<unsigned EIGHT_BYTE_INT>; |
|
|
72 template class octave_sort<vec_index *>; |
|
4850
|
73 #else |
|
|
74 struct vec_index |
|
|
75 { |
|
|
76 double vec; |
|
|
77 int indx; |
|
|
78 }; |
|
2928
|
79 |
|
4850
|
80 bool |
|
|
81 double_compare (double a, double b) |
|
|
82 { |
|
4991
|
83 return (xisnan (b) || (a < b)); |
|
4850
|
84 } |
|
2928
|
85 |
|
4850
|
86 bool |
|
|
87 double_compare (vec_index *a, vec_index *b) |
|
|
88 { |
|
4991
|
89 return (xisnan (b->vec) || (a->vec < b->vec)); |
|
4850
|
90 } |
|
|
91 |
|
4883
|
92 template class octave_sort<double>; |
|
|
93 template class octave_sort<vec_index *>; |
|
4850
|
94 #endif |
|
2928
|
95 |
|
4850
|
96 struct complex_vec_index |
|
|
97 { |
|
|
98 Complex vec; |
|
|
99 int indx; |
|
|
100 }; |
|
2928
|
101 |
|
4850
|
102 bool |
|
|
103 complex_compare (complex_vec_index *a, complex_vec_index *b) |
|
|
104 { |
|
4991
|
105 return (xisnan (b->vec) || (abs (a->vec) < abs (b->vec))); |
|
4850
|
106 } |
|
2928
|
107 |
|
4883
|
108 template class octave_sort<complex_vec_index *>; |
|
2928
|
109 |
|
|
110 static octave_value_list |
|
4850
|
111 mx_sort (NDArray &m, bool return_idx, int dim) |
|
2928
|
112 { |
|
|
113 octave_value_list retval; |
|
|
114 |
|
4850
|
115 if (m.length () < 1) |
|
|
116 return retval; |
|
|
117 |
|
|
118 dim_vector dv = m.dims (); |
|
4851
|
119 unsigned int ns = dv (dim); |
|
|
120 unsigned int iter = dv.numel () / ns; |
|
|
121 unsigned int stride = 1; |
|
|
122 for (unsigned int i = 0; i < (unsigned int)dim; i++) |
|
4850
|
123 stride *= dv(i); |
|
|
124 |
|
4853
|
125 #if defined (HAVE_IEEE754_DATA_FORMAT) && defined (EIGHT_BYTE_INT) |
|
4850
|
126 double *v = m.fortran_vec (); |
|
|
127 |
|
|
128 unsigned EIGHT_BYTE_INT *p = (unsigned EIGHT_BYTE_INT *)v; |
|
|
129 |
|
|
130 if (return_idx) |
|
|
131 { |
|
|
132 octave_sort<vec_index *> indexed_ieee754_sort (ieee754_compare); |
|
|
133 |
|
|
134 OCTAVE_LOCAL_BUFFER (vec_index *, vi, ns); |
|
|
135 OCTAVE_LOCAL_BUFFER (vec_index, vix, ns); |
|
|
136 |
|
4851
|
137 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
138 vi[i] = &vix[i]; |
|
|
139 |
|
|
140 NDArray idx (dv); |
|
|
141 |
|
4851
|
142 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
143 { |
|
4851
|
144 unsigned int offset = j; |
|
|
145 unsigned int offset2 = 0; |
|
4850
|
146 while (offset >= stride) |
|
|
147 { |
|
|
148 offset -= stride; |
|
|
149 offset2++; |
|
|
150 } |
|
|
151 offset += offset2 * stride * ns; |
|
2928
|
152 |
|
4853
|
153 // Flip the data in the vector so that int compares on |
|
|
154 // IEEE754 give the correct ordering. |
|
|
155 |
|
4851
|
156 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
157 { |
|
|
158 vi[i]->vec = FloatFlip (p[i*stride + offset]); |
|
|
159 vi[i]->indx = i + 1; |
|
|
160 } |
|
|
161 |
|
|
162 indexed_ieee754_sort.sort (vi, ns); |
|
|
163 |
|
4853
|
164 // Flip the data out of the vector so that int compares on |
|
|
165 // IEEE754 give the correct ordering |
|
|
166 |
|
4851
|
167 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
168 { |
|
|
169 p[i*stride + offset] = IFloatFlip (vi[i]->vec); |
|
|
170 idx(i*stride + offset) = vi[i]->indx; |
|
|
171 } |
|
2928
|
172 |
|
4853
|
173 // There are two representations of NaN. One will be sorted |
|
|
174 // to the beginning of the vector and the other to the end. |
|
|
175 // If it will be sorted to the beginning, fix things up. |
|
|
176 |
|
4850
|
177 if (lo_ieee_signbit (octave_NaN)) |
|
|
178 { |
|
|
179 unsigned int i = 0; |
|
|
180 while (xisnan(v[i++*stride+offset]) && i < ns); |
|
|
181 OCTAVE_LOCAL_BUFFER (double, itmp, i - 1); |
|
|
182 for (unsigned int l = 0; l < i -1; l++) |
|
|
183 itmp[l] = idx(l*stride + offset); |
|
|
184 for (unsigned int l = 0; l < ns - i + 1; l++) |
|
|
185 { |
|
|
186 v[l*stride + offset] = v[(l+i-1)*stride + offset]; |
|
|
187 idx(l*stride + offset) = idx((l+i-1)*stride + offset); |
|
|
188 } |
|
|
189 for (unsigned int k = 0, l = ns - i + 1; l < ns; l++, k++) |
|
|
190 { |
|
|
191 v[l*stride + offset] = octave_NaN; |
|
|
192 idx(l*stride + offset) = itmp[k]; |
|
|
193 } |
|
|
194 } |
|
|
195 } |
|
2928
|
196 |
|
4991
|
197 retval(1) = idx; |
|
2928
|
198 } |
|
4850
|
199 else |
|
2928
|
200 { |
|
4850
|
201 octave_sort<unsigned EIGHT_BYTE_INT> ieee754_sort; |
|
|
202 |
|
|
203 if (stride == 1) |
|
2928
|
204 { |
|
4851
|
205 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
206 { |
|
4853
|
207 // Flip the data in the vector so that int compares on |
|
|
208 // IEEE754 give the correct ordering. |
|
|
209 |
|
4851
|
210 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
211 p[i] = FloatFlip (p[i]); |
|
|
212 |
|
|
213 ieee754_sort.sort (p, ns); |
|
|
214 |
|
4853
|
215 // Flip the data out of the vector so that int compares |
|
|
216 // on IEEE754 give the correct ordering. |
|
|
217 |
|
4851
|
218 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
219 p[i] = IFloatFlip (p[i]); |
|
|
220 |
|
4853
|
221 // There are two representations of NaN. One will be |
|
|
222 // sorted to the beginning of the vector and the other |
|
|
223 // to the end. If it will be sorted to the beginning, |
|
|
224 // fix things up. |
|
|
225 |
|
4850
|
226 if (lo_ieee_signbit (octave_NaN)) |
|
|
227 { |
|
|
228 unsigned int i = 0; |
|
|
229 double *vtmp = (double *)p; |
|
|
230 while (xisnan(vtmp[i++]) && i < ns); |
|
|
231 for (unsigned int l = 0; l < ns - i + 1; l++) |
|
|
232 vtmp[l] = vtmp[l+i-1]; |
|
|
233 for (unsigned int l = ns - i + 1; l < ns; l++) |
|
|
234 vtmp[l] = octave_NaN; |
|
|
235 } |
|
|
236 |
|
|
237 p += ns; |
|
|
238 } |
|
2928
|
239 |
|
4850
|
240 } |
|
|
241 else |
|
|
242 { |
|
4851
|
243 OCTAVE_LOCAL_BUFFER (unsigned EIGHT_BYTE_INT, vi, ns); |
|
|
244 |
|
|
245 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
246 { |
|
4851
|
247 unsigned int offset = j; |
|
|
248 unsigned int offset2 = 0; |
|
4850
|
249 while (offset >= stride) |
|
|
250 { |
|
|
251 offset -= stride; |
|
|
252 offset2++; |
|
|
253 } |
|
|
254 offset += offset2 * stride * ns; |
|
|
255 |
|
4853
|
256 // Flip the data in the vector so that int compares on |
|
|
257 // IEEE754 give the correct ordering. |
|
|
258 |
|
4851
|
259 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
260 vi[i] = FloatFlip (p[i*stride + offset]); |
|
2928
|
261 |
|
4850
|
262 ieee754_sort.sort (vi, ns); |
|
|
263 |
|
4853
|
264 // Flip the data out of the vector so that int compares |
|
|
265 // on IEEE754 give the correct ordering. |
|
|
266 |
|
4851
|
267 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
268 p[i*stride + offset] = IFloatFlip (vi[i]); |
|
|
269 |
|
4853
|
270 // There are two representations of NaN. One will be |
|
|
271 // sorted to the beginning of the vector and the other |
|
|
272 // to the end. If it will be sorted to the beginning, |
|
|
273 // fix things up. |
|
|
274 |
|
4850
|
275 if (lo_ieee_signbit (octave_NaN)) |
|
|
276 { |
|
|
277 unsigned int i = 0; |
|
|
278 while (xisnan(v[i++*stride + offset]) && i < ns); |
|
|
279 for (unsigned int l = 0; l < ns - i + 1; l++) |
|
|
280 v[l*stride + offset] = v[(l+i-1)*stride + offset]; |
|
|
281 for (unsigned int l = ns - i + 1; l < ns; l++) |
|
|
282 v[l*stride + offset] = octave_NaN; |
|
|
283 } |
|
|
284 } |
|
2928
|
285 } |
|
|
286 } |
|
4850
|
287 #else |
|
|
288 if (return_idx) |
|
|
289 { |
|
|
290 double *v = m.fortran_vec (); |
|
|
291 octave_sort<vec_index *> indexed_double_sort (double_compare); |
|
2928
|
292 |
|
4850
|
293 OCTAVE_LOCAL_BUFFER (vec_index *, vi, ns); |
|
|
294 OCTAVE_LOCAL_BUFFER (vec_index, vix, ns); |
|
|
295 |
|
4851
|
296 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
297 vi[i] = &vix[i]; |
|
|
298 |
|
|
299 NDArray idx (dv); |
|
|
300 |
|
|
301 if (stride == 1) |
|
|
302 { |
|
4851
|
303 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
304 { |
|
4851
|
305 unsigned int offset = j * ns; |
|
4850
|
306 |
|
4851
|
307 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
308 { |
|
|
309 vi[i]->vec = v[i]; |
|
|
310 vi[i]->indx = i + 1; |
|
|
311 } |
|
|
312 |
|
|
313 indexed_double_sort.sort (vi, ns); |
|
|
314 |
|
4851
|
315 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
316 { |
|
|
317 v[i] = vi[i]->vec; |
|
|
318 idx(i + offset) = vi[i]->indx; |
|
|
319 } |
|
|
320 v += ns; |
|
|
321 } |
|
|
322 } |
|
|
323 else |
|
|
324 { |
|
4851
|
325 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
326 { |
|
4851
|
327 unsigned int offset = j; |
|
|
328 unsigned int offset2 = 0; |
|
4850
|
329 while (offset >= stride) |
|
|
330 { |
|
|
331 offset -= stride; |
|
|
332 offset2++; |
|
|
333 } |
|
|
334 offset += offset2 * stride * ns; |
|
|
335 |
|
4851
|
336 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
337 { |
|
|
338 vi[i]->vec = v[i*stride + offset]; |
|
|
339 vi[i]->indx = i + 1; |
|
|
340 } |
|
2928
|
341 |
|
4850
|
342 indexed_double_sort.sort (vi, ns); |
|
|
343 |
|
4851
|
344 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
345 { |
|
|
346 v[i*stride+offset] = vi[i]->vec; |
|
|
347 idx(i*stride+offset) = vi[i]->indx; |
|
|
348 } |
|
|
349 } |
|
|
350 } |
|
4991
|
351 retval(1) = idx; |
|
4850
|
352 } |
|
|
353 else |
|
|
354 { |
|
|
355 double *v = m.fortran_vec (); |
|
|
356 octave_sort<double> double_sort (double_compare); |
|
|
357 |
|
|
358 if (stride == 1) |
|
4851
|
359 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
360 { |
|
|
361 double_sort.sort (v, ns); |
|
|
362 v += ns; |
|
|
363 } |
|
|
364 else |
|
|
365 { |
|
|
366 OCTAVE_LOCAL_BUFFER (double, vi, ns); |
|
4851
|
367 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
368 { |
|
4851
|
369 unsigned int offset = j; |
|
|
370 unsigned int offset2 = 0; |
|
4850
|
371 while (offset >= stride) |
|
|
372 { |
|
|
373 offset -= stride; |
|
|
374 offset2++; |
|
|
375 } |
|
|
376 offset += offset2 * stride * ns; |
|
|
377 |
|
4851
|
378 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
379 vi[i] = v[i*stride + offset]; |
|
|
380 |
|
|
381 double_sort.sort (vi, ns); |
|
|
382 |
|
4851
|
383 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
384 v[i*stride + offset] = vi[i]; |
|
|
385 } |
|
|
386 } |
|
|
387 } |
|
|
388 #endif |
|
|
389 retval(0) = m; |
|
2928
|
390 return retval; |
|
|
391 } |
|
|
392 |
|
|
393 static octave_value_list |
|
4850
|
394 mx_sort (ComplexNDArray &m, bool return_idx, int dim) |
|
2928
|
395 { |
|
|
396 octave_value_list retval; |
|
|
397 |
|
4850
|
398 if (m.length () < 1) |
|
|
399 return retval; |
|
|
400 |
|
|
401 dim_vector dv = m.dims (); |
|
4851
|
402 unsigned int ns = dv (dim); |
|
|
403 unsigned int iter = dv.numel () / ns; |
|
|
404 unsigned int stride = 1; |
|
|
405 for (unsigned int i = 0; i < (unsigned int)dim; i++) |
|
4850
|
406 stride *= dv(i); |
|
2928
|
407 |
|
4850
|
408 octave_sort<complex_vec_index *> indexed_double_sort (complex_compare); |
|
|
409 |
|
|
410 Complex *v = m.fortran_vec (); |
|
2928
|
411 |
|
4850
|
412 OCTAVE_LOCAL_BUFFER (complex_vec_index *, vi, ns); |
|
|
413 OCTAVE_LOCAL_BUFFER (complex_vec_index, vix, ns); |
|
|
414 |
|
4851
|
415 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
416 vi[i] = &vix[i]; |
|
2928
|
417 |
|
4850
|
418 NDArray idx (dv); |
|
|
419 |
|
|
420 if (stride == 1) |
|
2928
|
421 { |
|
4851
|
422 for (unsigned int j = 0; j < iter; j++) |
|
2928
|
423 { |
|
4851
|
424 unsigned int offset = j * ns; |
|
2928
|
425 |
|
4851
|
426 for (unsigned int i = 0; i < ns; i++) |
|
4153
|
427 { |
|
4850
|
428 vi[i]->vec = v[i]; |
|
|
429 vi[i]->indx = i + 1; |
|
|
430 } |
|
|
431 |
|
|
432 indexed_double_sort.sort (vi, ns); |
|
|
433 |
|
|
434 if (return_idx) |
|
|
435 { |
|
4851
|
436 for (unsigned int i = 0; i < ns; i++) |
|
4153
|
437 { |
|
4850
|
438 v[i] = vi[i]->vec; |
|
|
439 idx(i + offset) = vi[i]->indx; |
|
4153
|
440 } |
|
|
441 } |
|
4850
|
442 else |
|
|
443 { |
|
4851
|
444 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
445 v[i] = vi[i]->vec; |
|
|
446 } |
|
|
447 v += ns; |
|
|
448 } |
|
|
449 } |
|
|
450 else |
|
|
451 { |
|
4851
|
452 for (unsigned int j = 0; j < iter; j++) |
|
4850
|
453 { |
|
4851
|
454 unsigned int offset = j; |
|
|
455 unsigned int offset2 = 0; |
|
4850
|
456 while (offset >= stride) |
|
|
457 { |
|
|
458 offset -= stride; |
|
|
459 offset2++; |
|
|
460 } |
|
|
461 offset += offset2 * stride * ns; |
|
2928
|
462 |
|
4851
|
463 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
464 { |
|
|
465 vi[i]->vec = v[i*stride + offset]; |
|
|
466 vi[i]->indx = i + 1; |
|
|
467 } |
|
|
468 |
|
|
469 indexed_double_sort.sort (vi, ns); |
|
|
470 |
|
|
471 if (return_idx) |
|
|
472 { |
|
4851
|
473 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
474 { |
|
|
475 v[i*stride + offset] = vi[i]->vec; |
|
|
476 idx(i*stride + offset) = vi[i]->indx; |
|
|
477 } |
|
|
478 } |
|
|
479 else |
|
|
480 { |
|
4851
|
481 for (unsigned int i = 0; i < ns; i++) |
|
4850
|
482 v[i*stride + offset] = vi[i]->vec; |
|
|
483 } |
|
2928
|
484 } |
|
|
485 } |
|
|
486 |
|
4850
|
487 if (return_idx) |
|
4991
|
488 retval(1) = idx; |
|
4850
|
489 |
|
|
490 retval(0) = m; |
|
2928
|
491 |
|
|
492 return retval; |
|
|
493 } |
|
|
494 |
|
4991
|
495 struct char_vec_index |
|
|
496 { |
|
|
497 char vec; |
|
|
498 int indx; |
|
|
499 }; |
|
|
500 |
|
|
501 bool |
|
|
502 char_compare (char_vec_index *a, char_vec_index *b) |
|
|
503 { |
|
|
504 return (a->vec < b->vec); |
|
|
505 } |
|
|
506 |
|
|
507 template class octave_sort<char>; |
|
|
508 template class octave_sort<char_vec_index *>; |
|
|
509 |
|
|
510 static octave_value_list |
|
|
511 mx_sort (charNDArray &m, bool return_idx, int dim) |
|
|
512 { |
|
|
513 octave_value_list retval; |
|
|
514 |
|
|
515 if (m.length () < 1) |
|
|
516 return retval; |
|
|
517 |
|
|
518 dim_vector dv = m.dims (); |
|
|
519 unsigned int ns = dv (dim); |
|
|
520 unsigned int iter = dv.numel () / ns; |
|
|
521 unsigned int stride = 1; |
|
|
522 for (unsigned int i = 0; i < static_cast <unsigned int> (dim); i++) |
|
|
523 stride *= dv(i); |
|
|
524 |
|
|
525 if (return_idx) |
|
|
526 { |
|
|
527 char *v = m.fortran_vec (); |
|
|
528 octave_sort<char_vec_index *> indexed_char_sort (char_compare); |
|
|
529 |
|
|
530 OCTAVE_LOCAL_BUFFER (char_vec_index *, vi, ns); |
|
|
531 OCTAVE_LOCAL_BUFFER (char_vec_index, vix, ns); |
|
|
532 |
|
|
533 for (unsigned int i = 0; i < ns; i++) |
|
|
534 vi[i] = &vix[i]; |
|
|
535 |
|
|
536 NDArray idx (dv); |
|
|
537 |
|
|
538 if (stride == 1) |
|
|
539 { |
|
|
540 for (unsigned int j = 0; j < iter; j++) |
|
|
541 { |
|
|
542 unsigned int offset = j * ns; |
|
|
543 |
|
|
544 for (unsigned int i = 0; i < ns; i++) |
|
|
545 { |
|
|
546 vi[i]->vec = v[i]; |
|
|
547 vi[i]->indx = i + 1; |
|
|
548 } |
|
|
549 |
|
|
550 indexed_char_sort.sort (vi, ns); |
|
|
551 |
|
|
552 for (unsigned int i = 0; i < ns; i++) |
|
|
553 { |
|
|
554 v[i] = vi[i]->vec; |
|
|
555 idx(i + offset) = vi[i]->indx; |
|
|
556 } |
|
|
557 v += ns; |
|
|
558 } |
|
|
559 } |
|
|
560 else |
|
|
561 { |
|
|
562 for (unsigned int j = 0; j < iter; j++) |
|
|
563 { |
|
|
564 unsigned int offset = j; |
|
|
565 unsigned int offset2 = 0; |
|
|
566 while (offset >= stride) |
|
|
567 { |
|
|
568 offset -= stride; |
|
|
569 offset2++; |
|
|
570 } |
|
|
571 offset += offset2 * stride * ns; |
|
|
572 |
|
|
573 for (unsigned int i = 0; i < ns; i++) |
|
|
574 { |
|
|
575 vi[i]->vec = v[i*stride + offset]; |
|
|
576 vi[i]->indx = i + 1; |
|
|
577 } |
|
|
578 |
|
|
579 indexed_char_sort.sort (vi, ns); |
|
|
580 |
|
|
581 for (unsigned int i = 0; i < ns; i++) |
|
|
582 { |
|
|
583 v[i*stride+offset] = vi[i]->vec; |
|
|
584 idx(i*stride+offset) = vi[i]->indx; |
|
|
585 } |
|
|
586 } |
|
|
587 } |
|
|
588 retval(1) = idx; |
|
|
589 } |
|
|
590 else |
|
|
591 { |
|
|
592 char *v = m.fortran_vec (); |
|
|
593 octave_sort<char> char_sort; |
|
|
594 |
|
|
595 if (stride == 1) |
|
|
596 for (unsigned int j = 0; j < iter; j++) |
|
|
597 { |
|
|
598 char_sort.sort (v, ns); |
|
|
599 v += ns; |
|
|
600 } |
|
|
601 else |
|
|
602 { |
|
|
603 OCTAVE_LOCAL_BUFFER (char, vi, ns); |
|
|
604 for (unsigned int j = 0; j < iter; j++) |
|
|
605 { |
|
|
606 unsigned int offset = j; |
|
|
607 unsigned int offset2 = 0; |
|
|
608 while (offset >= stride) |
|
|
609 { |
|
|
610 offset -= stride; |
|
|
611 offset2++; |
|
|
612 } |
|
|
613 offset += offset2 * stride * ns; |
|
|
614 |
|
|
615 for (unsigned int i = 0; i < ns; i++) |
|
|
616 vi[i] = v[i*stride + offset]; |
|
|
617 |
|
|
618 char_sort.sort (vi, ns); |
|
|
619 |
|
|
620 for (unsigned int i = 0; i < ns; i++) |
|
|
621 v[i*stride + offset] = vi[i]; |
|
|
622 } |
|
|
623 } |
|
|
624 } |
|
|
625 |
|
|
626 retval(0) = octave_value (m, true); |
|
|
627 return retval; |
|
|
628 } |
|
|
629 |
|
2928
|
630 DEFUN_DLD (sort, args, nargout, |
|
3369
|
631 "-*- texinfo -*-\n\ |
|
|
632 @deftypefn {Loadable Function} {[@var{s}, @var{i}] =} sort (@var{x})\n\ |
|
4850
|
633 @deftypefnx {Loadable Function} {[@var{s}, @var{i}] =} sort (@var{x}, @var{dim})\n\ |
|
3369
|
634 Return a copy of @var{x} with the elements elements arranged in\n\ |
|
|
635 increasing order. For matrices, @code{sort} orders the elements in each\n\ |
|
|
636 column.\n\ |
|
|
637 \n\ |
|
|
638 For example,\n\ |
|
|
639 \n\ |
|
|
640 @example\n\ |
|
|
641 @group\n\ |
|
|
642 sort ([1, 2; 2, 3; 3, 1])\n\ |
|
|
643 @result{} 1 1\n\ |
|
|
644 2 2\n\ |
|
|
645 3 3\n\ |
|
|
646 @end group\n\ |
|
|
647 @end example\n\ |
|
2928
|
648 \n\ |
|
3369
|
649 The @code{sort} function may also be used to produce a matrix\n\ |
|
|
650 containing the original row indices of the elements in the sorted\n\ |
|
|
651 matrix. For example,\n\ |
|
|
652 \n\ |
|
|
653 @example\n\ |
|
|
654 @group\n\ |
|
|
655 [s, i] = sort ([1, 2; 2, 3; 3, 1])\n\ |
|
|
656 @result{} s = 1 1\n\ |
|
|
657 2 2\n\ |
|
|
658 3 3\n\ |
|
|
659 @result{} i = 1 3\n\ |
|
|
660 2 1\n\ |
|
|
661 3 2\n\ |
|
|
662 @end group\n\ |
|
|
663 @end example\n\ |
|
4850
|
664 \n\ |
|
|
665 If the optional argument @var{dim} is given, then the matrix is sorted\n\ |
|
|
666 along the dimension defined by @var{dim}.\n\ |
|
|
667 \n\ |
|
|
668 For equal elements, the indices are such that the equal elements are listed\n\ |
|
|
669 in the order that appeared in the original list.\n\ |
|
|
670 \n\ |
|
|
671 The algorithm used in @code{sort} is optimized for the sorting of partially\n\ |
|
|
672 ordered lists.\n\ |
|
3369
|
673 @end deftypefn") |
|
2928
|
674 { |
|
|
675 octave_value_list retval; |
|
|
676 |
|
|
677 int nargin = args.length (); |
|
|
678 |
|
4850
|
679 if (nargin != 1 && nargin != 2) |
|
2928
|
680 { |
|
|
681 print_usage ("sort"); |
|
|
682 return retval; |
|
|
683 } |
|
|
684 |
|
4850
|
685 bool return_idx = nargout > 1; |
|
2928
|
686 |
|
|
687 octave_value arg = args(0); |
|
|
688 |
|
4850
|
689 int dim = 0; |
|
|
690 if (nargin == 2) |
|
|
691 dim = args(1).nint_value () - 1; |
|
|
692 |
|
|
693 dim_vector dv = ((const octave_complex_matrix&) arg) .dims (); |
|
|
694 if (error_state) |
|
|
695 { |
|
|
696 gripe_wrong_type_arg ("sort", arg); |
|
|
697 return retval; |
|
|
698 } |
|
|
699 if (nargin != 2) |
|
|
700 { |
|
|
701 // Find first non singleton dimension |
|
|
702 for (int i = 0; i < dv.length (); i++) |
|
|
703 if (dv(i) > 1) |
|
|
704 { |
|
|
705 dim = i; |
|
|
706 break; |
|
|
707 } |
|
|
708 } |
|
|
709 else |
|
|
710 { |
|
|
711 if (dim < 0 || dim > dv.length () - 1) |
|
|
712 { |
|
|
713 error ("sort: dim must be a valid dimension"); |
|
|
714 return retval; |
|
|
715 } |
|
|
716 } |
|
|
717 |
|
2928
|
718 if (arg.is_real_type ()) |
|
|
719 { |
|
4850
|
720 NDArray m = arg.array_value (); |
|
2928
|
721 |
|
|
722 if (! error_state) |
|
4850
|
723 retval = mx_sort (m, return_idx, dim); |
|
2928
|
724 } |
|
|
725 else if (arg.is_complex_type ()) |
|
|
726 { |
|
4850
|
727 ComplexNDArray cm = arg.complex_array_value (); |
|
2928
|
728 |
|
|
729 if (! error_state) |
|
4850
|
730 retval = mx_sort (cm, return_idx, dim); |
|
2928
|
731 } |
|
4991
|
732 else if (arg.is_string ()) |
|
|
733 { |
|
|
734 charNDArray chm = arg.char_array_value (); |
|
|
735 |
|
|
736 if (! error_state) |
|
|
737 retval = mx_sort (chm, return_idx, dim); |
|
|
738 } |
|
2928
|
739 else |
|
|
740 gripe_wrong_type_arg ("sort", arg); |
|
|
741 |
|
|
742 return retval; |
|
|
743 } |
|
|
744 |
|
|
745 /* |
|
|
746 ;;; Local Variables: *** |
|
|
747 ;;; mode: C++ *** |
|
|
748 ;;; End: *** |
|
|
749 */ |
