source: branches/hpc33/src/STIdxIter.cpp@ 2537

Last change on this file since 2537 was 2537, checked in by Takeshi Nakazato, 12 years ago

New Development: No

JIRA Issue: No

Ready for Test: Yes

Interface Changes: Yes

What Interface Changed: see below

Test Programs: List test programs

Put in Release Notes: Yes/No

Module(s): Module Names change impacts.

Description: Describe your changes here...

Defined STIdxIterExAcc class that is able to handle Int and String scalar
columns. Also, fixed a bug in ArrayIndexIterartorAcc.


File size: 11.8 KB
Line 
1#include <iostream>
2#include <casa/Utilities/GenSort.h>
3#include <casa/Arrays/Vector.h>
4#include <casa/Arrays/Matrix.h>
5#include <casa/Arrays/ArrayMath.h>
6#include <casa/Arrays/ArrayIO.h>
7#include <tables/Tables/ScalarColumn.h>
8#include "STIdxIter.h"
9
10namespace asap {
11
12IndexIterator::IndexIterator( IPosition &shape )
13 : niter_m( 0 )
14{
15 nfield_m = shape.nelements() ;
16 prod_m.resize( nfield_m ) ;
17 idx_m.resize( nfield_m ) ;
18 prod_m[0] = 1 ;
19 idx_m[0] = 0 ;
20 for ( uInt i = 1 ; i < nfield_m ; i++ ) {
21 prod_m[i] = shape[i-1] * prod_m[i-1] ;
22 idx_m[i] = 0 ;
23 }
24 maxiter_m = prod_m[nfield_m-1] * shape[nfield_m-1] ;
25// cout << "maxiter_m=" << maxiter_m << endl ;
26}
27
28Bool IndexIterator::pastEnd()
29{
30 return niter_m >= maxiter_m ;
31}
32
33void IndexIterator::next()
34{
35 niter_m++ ;
36 uInt x = niter_m ;
37 for ( Int i = nfield_m-1 ; i >= 0 ; i-- ) {
38 idx_m[i] = x / prod_m[i] ;
39 //cout << "i=" << i << ": prod=" << prod_m[i]
40 // << ", idx=" << idx_m[i] << endl ;
41 x %= prod_m[i] ;
42 }
43}
44
45ArrayIndexIterator::ArrayIndexIterator( Matrix<uInt> &arr,
46 vector< vector<uInt> > idlist )
47 : arr_m( arr ),
48 pos_m( 1 )
49{
50 ncol_m = arr_m.ncolumn() ;
51 nrow_m = arr_m.nrow() ;
52 vector< vector<uInt> > l = idlist ;
53 if ( l.size() != ncol_m ) {
54 l.resize( ncol_m ) ;
55 for ( uInt i = 0 ; i < ncol_m ; i++ ) {
56 Vector<uInt> a( arr_m.column( i ).copy() ) ;
57 uInt n = genSort( a, Sort::Ascending, Sort::QuickSort|Sort::NoDuplicates ) ;
58 a.resize(n,True) ;
59 a.tovector( l[i] ) ;
60 }
61 }
62 idxlist_m = l ;
63 IPosition shape( ncol_m ) ;
64 for ( uInt i = 0 ; i < ncol_m ; i++ ) {
65 shape[i] = idxlist_m[i].size() ;
66 }
67// cout << "shape=" << shape << endl ;
68 iter_m = new IndexIterator( shape ) ;
69 current_m.resize( ncol_m ) ;
70}
71
72ArrayIndexIterator::~ArrayIndexIterator()
73{
74 delete iter_m ;
75}
76
77Vector<uInt> ArrayIndexIterator::current()
78{
79 Block<uInt> idx = iter_m->current() ;
80 for ( uInt i = 0 ; i < ncol_m ; i++ )
81 current_m[i] = idxlist_m[i][idx[i]] ;
82 return current_m ;
83}
84
85Bool ArrayIndexIterator::pastEnd()
86{
87 return iter_m->pastEnd() ;
88}
89
90ArrayIndexIteratorNormal::ArrayIndexIteratorNormal( Matrix<uInt> &arr,
91 vector< vector<uInt> > idlist )
92 : ArrayIndexIterator( arr, idlist )
93{
94 storage_m.resize( nrow_m ) ;
95}
96
97void ArrayIndexIteratorNormal::next()
98{
99 iter_m->next() ;
100}
101
102Vector<uInt> ArrayIndexIteratorNormal::getRows( StorageInitPolicy policy )
103{
104 Vector<uInt> v = current() ;
105 uInt len = 0 ;
106 uInt *p = storage_m.storage() ;
107 for ( uInt i = 0 ; i < nrow_m ; i++ ) {
108 if ( allEQ( v, arr_m.row( i ) ) ) {
109 *p = i ;
110 len++ ;
111 p++ ;
112 }
113 }
114 pos_m[0] = len ;
115 p = storage_m.storage() ;
116 return Vector<uInt>( pos_m, p, policy ) ;
117}
118
119ArrayIndexIteratorAcc::ArrayIndexIteratorAcc( Matrix<uInt> &arr,
120 vector< vector<uInt> > idlist )
121 : ArrayIndexIterator( arr, idlist )
122{
123 storage_m.resize( arr_m.nelements()+nrow_m ) ;
124 uInt *p = storage_m.storage() + arr_m.nelements() ;
125 for ( uInt i = 0 ; i < nrow_m ; i++ ) {
126 *p = i ;
127 p++ ;
128 }
129 len_m.resize( ncol_m+1 ) ;
130 p = len_m.storage() ;
131 for ( uInt i = 0 ; i < ncol_m+1 ; i++ ) {
132 *p = nrow_m ;
133 p++ ;
134// cout << "len[" << i << "]=" << len_m[i] << endl ;
135 }
136 prev_m = iter_m->current() ;
137 p = prev_m.storage() ;
138 for ( uInt i = 0 ; i < ncol_m ; i++ ) {
139 *p = *p - 1 ;
140 p++ ;
141 }
142// cout << "prev_m=" << Vector<uInt>(IPosition(1,ncol_m),prev_m.storage()) << endl ;
143 skip_m.resize( ncol_m ) ;
144 for ( uInt i = 0 ; i < ncol_m ; i++ ) {
145 skip_m[i] = (Bool)(idxlist_m[i].size()==1) ;
146// cout << "skip_m[" << i << "]=" << skip_m[i] << endl ;
147 }
148}
149
150void ArrayIndexIteratorAcc::next()
151{
152 prev_m = iter_m->current() ;
153 iter_m->next() ;
154}
155
156Vector<uInt> ArrayIndexIteratorAcc::getRows( StorageInitPolicy policy )
157{
158 Block<uInt> v = iter_m->current() ;
159 Int c = isChanged( v ) ;
160// cout << "c=" << c << endl ;
161 if ( c < 0 ) {
162 pos_m[0] = len_m[0] ;
163 uInt offset = 0 ;
164 while( offset < ncol_m && skip_m[offset] )
165 offset++ ;
166// cout << "offset=" << offset << endl ;
167 return Vector<uInt>( pos_m, storage_m.storage()+offset*nrow_m, policy ) ;
168 }
169 Int offset = c + 1 ;
170 while( offset < ncol_m && skip_m[offset] )
171 offset++ ;
172// cout << "offset = " << offset << endl ;
173 uInt *base = storage_m.storage() + offset * nrow_m ;
174// cout << "len_m[c+1]=" << len_m[c+1] << endl ;
175// cout << "base=" << Vector<uInt>(IPosition(1,len_m[c+1]),base,SHARE) << endl ;
176 for ( Int i = c ; i >= 0 ; i-- ) {
177 base = updateStorage( i, base, idxlist_m[i][v[i]] ) ;
178// cout << "len_m[" << i << "]=" << len_m[i] << endl ;
179// cout << "base=" << Vector<uInt>(IPosition(1,len_m[i]),base,SHARE) << endl ;
180 }
181 pos_m[0] = len_m[0] ;
182// cout << "pos_m=" << pos_m << endl ;
183 return Vector<uInt>( pos_m, base, policy ) ;
184}
185
186Int ArrayIndexIteratorAcc::isChanged( Block<uInt> &idx )
187{
188 Int i = ncol_m - 1 ;
189 while( i >= 0 && idx[i] == prev_m[i] )
190 i-- ;
191 return i ;
192}
193
194uInt *ArrayIndexIteratorAcc::updateStorage( Int &icol,
195 uInt *base,
196 uInt &v )
197{
198 uInt *p ;
199 if ( skip_m[icol] ) {
200 // skip update, just update len_m[icol] and pass appropriate pointer
201// cout << "skip " << icol << endl ;
202 p = base ;
203 len_m[icol] = len_m[icol+1] ;
204 }
205 else {
206// cout << "update " << icol << endl ;
207 uInt len = 0 ;
208 p = storage_m.storage() + icol * nrow_m ;
209 uInt *work = p ;
210 Bool b ;
211 const uInt *arr_p = arr_m.getStorage( b ) ;
212 long offset = 0 ;
213 // warr_p points a first element of (icol)-th column
214 const uInt *warr_p = arr_p + icol * nrow_m ;
215 for ( uInt i = 0 ; i < len_m[icol+1] ; i++ ) {
216 // increment warr_p by (*(base)-*(base-1))
217 warr_p += *base - offset ;
218 // check if target element is equal to value specified
219 if ( *warr_p == v ) {
220 // then, add current index to storage_m
221 // cout << "add " << *base << endl ;
222 *work = *base ;
223 len++ ;
224 work++ ;
225 }
226 // update offset
227 offset = *base ;
228 // next index
229 base++ ;
230 }
231 arr_m.freeStorage( arr_p, b ) ;
232 len_m[icol] = len ;
233 }
234 return p ;
235}
236
237STIdxIter::STIdxIter()
238{
239 iter_m = 0 ;
240}
241
242STIdxIter::STIdxIter( const string &name,
243 const vector<string> &cols )
244{
245}
246
247STIdxIter::STIdxIter( const CountedPtr<Scantable> &s,
248 const vector<string> &cols )
249{
250}
251
252STIdxIter::~STIdxIter()
253{
254 if ( iter_m != 0 )
255 delete iter_m ;
256}
257
258vector<uInt> STIdxIter::tovector( Vector<uInt> v )
259{
260 vector<uInt> ret ;
261 v.tovector( ret ) ;
262 return ret ;
263}
264
265Vector<uInt> STIdxIter::getRows( StorageInitPolicy policy )
266{
267 return iter_m->getRows( policy ) ;
268}
269
270STIdxIterNormal::STIdxIterNormal()
271 : STIdxIter()
272{
273}
274
275STIdxIterNormal::STIdxIterNormal( const string &name,
276 const vector<string> &cols )
277 : STIdxIter( name, cols )
278{
279 Table t( name, Table::Old ) ;
280 init( t, cols ) ;
281}
282
283STIdxIterNormal::STIdxIterNormal( const CountedPtr<Scantable> &s,
284 const vector<string> &cols )
285 : STIdxIter( s, cols )
286{
287 init( s->table(), cols ) ;
288}
289
290STIdxIterNormal::~STIdxIterNormal()
291{
292}
293
294void STIdxIterNormal::init( Table &t,
295 const vector<string> &cols )
296{
297 uInt ncol = cols.size() ;
298 uInt nrow = t.nrow() ;
299 Matrix<uInt> arr( nrow, ncol ) ;
300 ROScalarColumn<uInt> col ;
301 Vector<uInt> v ;
302 for ( uInt i = 0 ; i < ncol ; i++ ) {
303 col.attach( t, cols[i] ) ;
304 v.reference( arr.column( i ) ) ;
305 col.getColumn( v ) ;
306 }
307 iter_m = new ArrayIndexIteratorNormal( arr ) ;
308}
309
310STIdxIterAcc::STIdxIterAcc()
311 : STIdxIter()
312{
313}
314
315STIdxIterAcc::STIdxIterAcc( const string &name,
316 const vector<string> &cols )
317 : STIdxIter( name, cols )
318{
319 Table t( name, Table::Old ) ;
320 init( t, cols ) ;
321}
322
323STIdxIterAcc::STIdxIterAcc( const CountedPtr<Scantable> &s,
324 const vector<string> &cols )
325 : STIdxIter( s, cols )
326{
327 init( s->table(), cols ) ;
328}
329
330STIdxIterAcc::~STIdxIterAcc()
331{
332}
333
334void STIdxIterAcc::init( Table &t,
335 const vector<string> &cols )
336{
337 uInt ncol = cols.size() ;
338 uInt nrow = t.nrow() ;
339 // array shape here is as follows if cols=["BEAMNO","POLNO","IFNO"]:
340 // [[B0,B1,B2,...,BN],
341 // [P0,P1,P2,...,PN],
342 // [I0,I1,I2,...,IN]]
343 // order of internal storage is
344 // [B0,B1,B2,..,BN,P0,P1,P2,...,PN,I0,I1,I2,...,IN]
345 Matrix<uInt> arr( nrow, ncol ) ;
346 Vector<uInt> v ;
347 ROScalarColumn<uInt> col ;
348 for ( uInt i = 0 ; i < ncol ; i++ ) {
349 col.attach( t, cols[i] ) ;
350 v.reference( arr.column( i ) ) ;
351 col.getColumn( v ) ;
352 }
353 iter_m = new ArrayIndexIteratorAcc( arr ) ;
354}
355
356STIdxIterExAcc::STIdxIterExAcc()
357 : STIdxIter()
358{
359}
360
361STIdxIterExAcc::STIdxIterExAcc( const string &name,
362 const vector<string> &cols )
363 : STIdxIter( name, cols )
364{
365 Table t( name, Table::Old ) ;
366 init( t, cols ) ;
367}
368
369STIdxIterExAcc::STIdxIterExAcc( const CountedPtr<Scantable> &s,
370 const vector<string> &cols )
371 : STIdxIter( s, cols )
372{
373 init( s->table(), cols ) ;
374}
375
376STIdxIterExAcc::~STIdxIterExAcc()
377{
378}
379
380void STIdxIterExAcc::init( Table &t,
381 const vector<string> &cols )
382{
383 uInt ncol = cols.size() ;
384 uInt nrow = t.nrow() ;
385 // array shape here is as follows if cols=["BEAMNO","POLNO","IFNO"]:
386 // [[B0,B1,B2,...,BN],
387 // [P0,P1,P2,...,PN],
388 // [I0,I1,I2,...,IN]]
389 // order of internal storage is
390 // [B0,B1,B2,..,BN,P0,P1,P2,...,PN,I0,I1,I2,...,IN]
391 Matrix<uInt> arr( nrow, ncol ) ;
392 Vector<uInt> v ;
393 ROScalarColumn<uInt> col ;
394 ROScalarColumn<String> strCol ;
395 ROScalarColumn<Int> intCol ;
396 for ( uInt i = 0 ; i < ncol ; i++ ) {
397 v.reference( arr.column( i ) ) ;
398 if ( cols[i] == "SRCTYPE" ) {
399 intCol.attach( t, cols[i] ) ;
400 Vector<Int> srctype = intCol.getColumn() ;
401 processIntCol( srctype, v, srctype_m ) ;
402 srctypeid_m = i ;
403 }
404 else if ( cols[i] == "SRCNAME" ) {
405 strCol.attach( t, cols[i] ) ;
406 Vector<String> srcname = strCol.getColumn() ;
407 processStrCol( srcname, v, srcname_m ) ;
408 srcnameid_m = i ;
409 }
410 else {
411 col.attach( t, cols[i] ) ;
412 col.getColumn( v ) ;
413 }
414 }
415 iter_m = new ArrayIndexIteratorAcc( arr ) ;
416}
417
418void STIdxIterExAcc::processIntCol( Vector<Int> &in,
419 Vector<uInt> &out,
420 Block<Int> &val )
421{
422 uInt len = in.nelements() ;
423 Vector<Int> tmp = in.copy() ;
424 uInt n = genSort( tmp, Sort::Ascending, Sort::QuickSort|Sort::NoDuplicates ) ;
425 val.resize( n ) ;
426 for ( uInt i = 0 ; i < n ; i++ ) {
427 val[i] = tmp[i] ;
428// cout << "val[" << i << "]=" << val[i] << endl ;
429 }
430 out.resize( len ) ;
431 for ( uInt i = 0 ; i < len ; i++ ) {
432 for ( uInt j = 0 ; j < n ; j++ ) {
433 if ( in[i] == val[j] ) {
434 out[i] = j ;
435 break ;
436 }
437 }
438 }
439}
440
441void STIdxIterExAcc::processStrCol( Vector<String> &in,
442 Vector<uInt> &out,
443 Block<String> &val )
444{
445 uInt len = in.nelements() ;
446 Vector<String> tmp = in.copy() ;
447 uInt n = genSort( tmp, Sort::Ascending, Sort::QuickSort|Sort::NoDuplicates ) ;
448 val.resize( n ) ;
449 for ( uInt i = 0 ; i < n ; i++ ) {
450 val[i] = tmp[i] ;
451// cout << "val[" << i << "]=" << val[i] << endl ;
452 }
453 out.resize( len ) ;
454 for ( uInt i = 0 ; i < len ; i++ ) {
455 for ( uInt j = 0 ; j < n ; j++ ) {
456 if ( in[i] == val[j] ) {
457 out[i] = j ;
458 break ;
459 }
460 }
461 }
462}
463
464Int STIdxIterExAcc::getSrcType()
465{
466 if ( srctype_m.nelements() > 0 )
467 return srctype_m[iter_m->current()[srctypeid_m]] ;
468 else
469 return -999 ;
470}
471
472String STIdxIterExAcc::getSrcName()
473{
474 if ( srcname_m.nelements() > 0 )
475 return srcname_m[iter_m->current()[srcnameid_m]] ;
476 else
477 return "" ;
478}
479
480} // namespace
Note: See TracBrowser for help on using the repository browser.