Changeset 2360
- Timestamp:
- 12/06/11 13:13:24 (13 years ago)
- Location:
- trunk
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/python/asapgrid.py
r2358 r2360 13 13 self.gridder._setin( infile ) 14 14 15 def setPolList( self, pollist ): 16 self.gridder._setpollist( pollist ) 17 15 18 def defineImage( self, nx=-1, ny=-1, cellx='', celly='', center='' ): 16 19 self.gridder._defineimage( nx, ny, cellx, celly, center ) … … 25 28 self.outfile = self.gridder._save( outfile ) 26 29 27 def plot( self, plotchan=-1 ):30 def plot( self, plotchan=-1, plotpol=-1 ): 28 31 plotter = _SDGridPlotter( self.infile, self.outfile ) 29 plotter.plot( chan=plotchan )32 plotter.plot( chan=plotchan, pol=plotpol ) 30 33 31 34 class _SDGridPlotter: … … 41 44 self.ny = -1 42 45 self.nchan = 0 46 self.npol = 0 47 self.pollist = [] 43 48 self.cellx = 0.0 44 49 self.celly = 0.0 … … 57 62 del s 58 63 59 idx = spectra.nonzero()[1]60 #self.nonzero = self.pointing.take( idx, axis=1 )61 62 64 s = scantable( self.outfile, average=False ) 63 self.grid = numpy.array( s.get_directionval() ).transpose()64 dirstring = s.get_direction()65 65 nrow = s.nrow() 66 spectra = []66 pols = numpy.ones( nrow, dtype=int ) 67 67 for i in xrange(nrow): 68 spectra.append( s._getspectrum( i ) ) 69 spectra = numpy.array( spectra ).transpose() 68 pols[i] = s.getpol(i) 69 self.pollist, indices = numpy.unique( pols, return_inverse=True ) 70 self.npol = len(self.pollist) 71 self.pollist = self.pollist[indices[:self.npol]] 72 #print 'pollist=',self.pollist 73 #print 'npol=',self.npol 74 #print 'nrow=',nrow 75 dirstring = numpy.array(s.get_direction()).take(range(0,nrow,self.npol)) 76 self.grid = numpy.array( s.get_directionval() ).take(range(0,nrow,self.npol),axis=0).transpose() 77 spectra = numpy.zeros( (self.npol,self.nchan,nrow/self.npol), dtype=float ) 78 irow = 0 79 for i in xrange(nrow/self.npol): 80 for ip in xrange(self.npol): 81 spectra[ip,:,i] = s._getspectrum( irow ) 82 irow += 1 70 83 71 84 idx = 0 … … 73 86 while ( dirstring[idx].split()[-1] == d0 ): 74 87 idx += 1 75 76 88 self.ny = idx 77 self.nx = nrow / idx 89 self.nx = nrow / (self.npol * idx ) 90 #print 'nx,ny=',self.nx,self.ny 78 91 79 92 self.cellx = abs( self.grid[0][0] - self.grid[0][1] ) 80 93 self.celly = abs( self.grid[1][0] - self.grid[1][self.ny] ) 94 #print 'cellx,celly=',self.cellx,self.celly 81 95 82 self.data = spectra.reshape( (self.n chan,self.nx,self.ny) )96 self.data = spectra.reshape( (self.npol,self.nchan,self.nx,self.ny) ) 83 97 84 def plot( self, chan=-1 ): 98 def plot( self, chan=-1, pol=-1 ): 99 if pol < 0: 100 data = self.data.mean(axis=0) 101 opt = 'averaged over pol' 102 else: 103 idx = self.pollist.tolist().index( pol ) 104 #print 'idx=',idx 105 data = self.data[idx] 106 opt = 'pol %s'%(pol) 85 107 if chan < 0: 86 data = self.data.mean(axis=0)87 title = 'Gridded Image (averaged over channel)'108 data = data.mean(axis=0) 109 opt += ', averaged over channel' 88 110 else: 89 data = self.data[chan] 90 title = 'Gridded Image (channel %s)'%(chan) 111 data = data[chan] 112 opt += ', channel %s'%(chan) 113 title = 'Gridded Image (%s)'%(opt) 91 114 pl.figure(10) 92 115 pl.clf() 93 116 pl.plot(self.grid[0],self.grid[1],'.',color='blue') 94 117 pl.plot(self.pointing[0],self.pointing[1],'.',color='red') 95 #pl.plot(self.nonzero[0],self.nonzero[1],'o',color='green')96 118 extent=[self.grid[0].min()-0.5*self.cellx, 97 119 self.grid[0].max()+0.5*self.cellx, -
trunk/src/STGrid.cpp
r2359 r2360 13 13 14 14 #include <tables/Tables/Table.h> 15 #include <tables/Tables/TableRecord.h> 16 #include <tables/Tables/ExprNode.h> 15 17 #include <tables/Tables/ScalarColumn.h> 16 18 #include <tables/Tables/ArrayColumn.h> … … 60 62 tab_ = Table( infile_ ) ; 61 63 } 64 } 65 66 void STGrid::setPolList( vector<unsigned int> pols ) 67 { 68 //pollist_ = Vector<uInt>( pols ) ; 69 pollist_.assign( Vector<uInt>( pols ) ) ; 70 cout << "pollist_ = " << pollist_ << endl ; 62 71 } 63 72 … … 126 135 127 136 // retrieve data 128 Matrix<Float> spectra ;137 Cube<Float> spectra ; 129 138 Matrix<Double> direction ; 130 Matrix<uChar> flagtra ;131 Vector<uInt> rflag ;139 Cube<uChar> flagtra ; 140 Matrix<uInt> rflag ; 132 141 getData( infile_, spectra, direction, flagtra, rflag ) ; 133 142 IPosition sshape = spectra.shape() ; 134 Int nchan = sshape[0] ; 135 Int nrow = sshape[1] ; 136 //cout << "data.shape()=" << data.shape() << endl ; 143 Int npol = sshape[0] ; 144 Int nchan = sshape[1] ; 145 Int nrow = sshape[2] ; 146 cout << "spectra.shape()=" << spectra.shape() << endl ; 147 cout << "max(spectra) = " << max(spectra) << endl ; 137 148 138 149 // flagtra: uChar -> Int 139 150 // rflag: uInt -> Int 140 Matrix<Int> flagI ; 141 Vector<Int> rflagI ; 142 toInt( flagtra, flagI ) ; 143 toInt( rflag, rflagI ) ; 151 // Matrix<Int> flagI ; 152 // Vector<Int> rflagI ; 153 Cube<Int> flagI ; 154 Matrix<Int> rflagI ; 155 toInt( &flagtra, &flagI ) ; 156 toInt( &rflag, &rflagI ) ; 144 157 145 //cout << "flagI.shape() = " << flagI.shape() << endl ;146 //cout << "rflagI.shape() = " << rflagI.shape() << endl ;158 cout << "flagI.shape() = " << flagI.shape() << endl ; 159 cout << "rflagI.shape() = " << rflagI.shape() << endl ; 147 160 148 161 // grid parameter … … 177 190 Matrix<Float> weight( IPosition( 2, nchan, nrow ), 1.0 ) ; 178 191 179 // call ggridsd192 // // call ggridsd 180 193 Bool deletePos, deleteData, deleteWgt, deleteFlag, deleteFlagR, deleteConv, deleteDataG, deleteWgtG ; 181 194 Double *xypos_p = xypos.getStorage( deletePos ) ; 182 Matrix<Complex> dataC( spectra.shape(), 0.0 ) ; 195 //Matrix<Complex> dataC( spectra.shape(), 0.0 ) ; 196 Cube<Complex> dataC( spectra.shape(), 0.0 ) ; 183 197 setReal( dataC, spectra ) ; 184 198 const Complex *data_p = dataC.getStorage( deleteData ) ; … … 187 201 const Int *rflag_p = rflagI.getStorage( deleteFlagR ) ; 188 202 Float *conv_p = convFunc.getStorage( deleteConv ) ; 189 Int npol = 1 ;203 //Int npol = 1 ; 190 204 // Extend grid plane with convSupport_ 191 205 //IPosition gshape( 4, nx_, ny_, npol, nchan ) ; … … 235 249 gdata_p, 236 250 wdata_p, 237 //&nx_,238 //&ny_,239 251 &gnx, 240 252 &gny, … … 381 393 382 394 void STGrid::getData( String &infile, 383 Matrix<Float> &spectra,395 Cube<Float> &spectra, 384 396 Matrix<Double> &direction, 385 Matrix<uChar> &flagtra,386 Vector<uInt> &rflag )397 Cube<uChar> &flagtra, 398 Matrix<uInt> &rflag ) 387 399 { 388 400 Table tab( infile ) ; 389 ROArrayColumn<Float> spectraCol( tab, "SPECTRA" ) ; 390 ROArrayColumn<Double> directionCol( tab, "DIRECTION" ) ; 391 ROArrayColumn<uChar> flagtraCol( tab, "FLAGTRA" ) ; 392 ROScalarColumn<uInt> rflagCol( tab, "FLAGROW" ) ; 393 spectraCol.getColumn( spectra ) ; 394 directionCol.getColumn( direction ) ; 395 flagtraCol.getColumn( flagtra ) ; 396 rflagCol.getColumn( rflag ) ; 397 } 398 399 void STGrid::toInt( Matrix<uChar> &u, Matrix<Int> &v ) 400 { 401 uInt len = u.nelements() ; 401 //uInt npol = tab.keywordSet().asuInt( "nPol" ) ; 402 ROScalarColumn<uInt> polnoCol( tab, "POLNO" ) ; 403 //uInt npol = max( polnoCol.getColumn() ) + 1 ; 404 Vector<uInt> pols = polnoCol.getColumn() ; 405 Vector<uInt> pollistOrg ; 406 uInt npolOrg = 0 ; 407 for ( uInt i = 0 ; i < pols.size() ; i++ ) { 408 if ( allNE( pollistOrg, pols[i] ) ) { 409 pollistOrg.resize( npolOrg+1, True ) ; 410 pollistOrg[npolOrg] = pols[i] ; 411 npolOrg++ ; 412 } 413 } 414 if ( pollist_.size() == 0 ) 415 pollist_ = pollistOrg ; 416 else { 417 Vector<uInt> newlist ; 418 uInt newsize = 0 ; 419 for ( uInt i = 0 ; i < pollist_.size() ; i++ ) { 420 if ( anyEQ( pols, pollist_[i] ) ) { 421 newlist.resize( newsize+1, True ) ; 422 newlist[newsize] = pollist_[i] ; 423 newsize++ ; 424 } 425 } 426 pollist_ = newlist ; 427 } 428 uInt npol = pollist_.size() ; 429 ROArrayColumn<uChar> tmpCol( tab, "FLAGTRA" ) ; 430 uInt nchan = tmpCol( 0 ).nelements() ; 431 uInt nrow = tab.nrow() / npolOrg ; 432 cout << "npol = " << npol << endl ; 433 cout << "nchan = " << nchan << endl ; 434 cout << "nrow = " << nrow << endl ; 435 spectra.resize( npol, nchan, nrow ) ; 436 flagtra.resize( npol, nchan, nrow ) ; 437 rflag.resize( npol, nrow ) ; 438 cout << "pollist_=" << pollist_ << endl ; 439 for ( uInt ipol = 0 ; ipol < npol ; ipol++ ) { 440 Table subt = tab( tab.col("POLNO") == pollist_[ipol] ) ; 441 ROArrayColumn<Float> spectraCol( subt, "SPECTRA" ) ; 442 ROArrayColumn<Double> directionCol( subt, "DIRECTION" ) ; 443 ROArrayColumn<uChar> flagtraCol( subt, "FLAGTRA" ) ; 444 ROScalarColumn<uInt> rflagCol( subt, "FLAGROW" ) ; 445 Matrix<Float> tmpF = spectra.yzPlane( ipol ) ; 446 Matrix<uChar> tmpUC = flagtra.yzPlane( ipol ) ; 447 Vector<uInt> tmpUI = rflag.row( ipol ) ; 448 cout << "tmpF.shape()=" << tmpF.shape() << endl ; 449 cout << "tmpUC.shape()=" << tmpUC.shape() << endl ; 450 cout << "tmpUI.shape()=" << tmpUI.shape() << endl ; 451 spectraCol.getColumn( tmpF ) ; 452 flagtraCol.getColumn( tmpUC ) ; 453 rflagCol.getColumn( tmpUI ) ; 454 if ( ipol == 0 ) 455 directionCol.getColumn( direction ) ; 456 } 457 cout << "getData end" << endl ; 458 } 459 460 void STGrid::toInt( Array<uChar> *u, Array<Int> *v ) 461 { 462 uInt len = u->nelements() ; 402 463 Int *int_p = new Int[len] ; 403 464 Bool deleteIt ; 404 const uChar *data_p = u .getStorage( deleteIt ) ;465 const uChar *data_p = u->getStorage( deleteIt ) ; 405 466 Int *i_p = int_p ; 406 467 const uChar *u_p = data_p ; … … 410 471 u_p++ ; 411 472 } 412 u .freeStorage( data_p, deleteIt ) ;413 v .takeStorage( u.shape(), int_p, TAKE_OVER ) ;414 } 415 416 void STGrid::toInt( Vector<uInt> &u, Vector<Int> &v )417 { 418 uInt len = u .nelements() ;473 u->freeStorage( data_p, deleteIt ) ; 474 v->takeStorage( u->shape(), int_p, TAKE_OVER ) ; 475 } 476 477 void STGrid::toInt( Array<uInt> *u, Array<Int> *v ) 478 { 479 uInt len = u->nelements() ; 419 480 Int *int_p = new Int[len] ; 420 481 Bool deleteIt ; 421 const uInt *data_p = u .getStorage( deleteIt ) ;482 const uInt *data_p = u->getStorage( deleteIt ) ; 422 483 Int *i_p = int_p ; 423 484 const uInt *u_p = data_p ; … … 427 488 u_p++ ; 428 489 } 429 u .freeStorage( data_p, deleteIt ) ;430 v .takeStorage( u.shape(), int_p, TAKE_OVER ) ;490 u->freeStorage( data_p, deleteIt ) ; 491 v->takeStorage( u->shape(), int_p, TAKE_OVER ) ; 431 492 } 432 493 … … 532 593 string STGrid::saveData( string outfile ) 533 594 { 534 Int polno = 0 ;595 //Int polno = 0 ; 535 596 string outfile_ ; 536 597 if ( outfile.size() == 0 ) { … … 550 611 CountedPtr<Scantable> out( new Scantable( *ref, True ) ) ; 551 612 Table tab = out->table() ; 552 Int nrow = nx_ * ny_ ;553 tab.addRow( nrow ) ;554 613 IPosition dshape = data_.shape() ; 555 614 Int npol = dshape[2] ; 556 615 Int nchan = dshape[3] ; 616 Int nrow = nx_ * ny_ * npol ; 617 tab.rwKeywordSet().define( "nPol", npol ) ; 618 tab.addRow( nrow ) ; 557 619 Vector<Double> cpix( 2 ) ; 558 620 cpix(0) = Double( nx_ - 1 ) * 0.5 ; … … 574 636 spectraCol.put( irow, sp ) ; 575 637 directionCol.put( irow, dir ) ; 576 polnoCol.put( irow, pol no) ;638 polnoCol.put( irow, pollist_[ipol] ) ; 577 639 irow++ ; 578 640 } -
trunk/src/STGrid.h
r2356 r2360 16 16 #include <fstream> 17 17 #include <string> 18 #include <vector> 18 19 19 20 #include <casa/BasicSL/String.h> 20 21 #include <casa/Arrays/Vector.h> 21 //#include <casa/Arrays/Matrix.h>22 //#include <casa/Arrays/Cube.h>22 #include <casa/Arrays/Matrix.h> 23 #include <casa/Arrays/Cube.h> 23 24 // #include <casa/Arrays/ArrayMath.h> 24 25 // #include <casa/Inputs/Input.h> … … 47 48 48 49 void setFileIn( const string infile ) ; 49 // void setFileOut( const string outfile ) ; 50 51 void setPolList( vector<unsigned int> pols ) ; 50 52 51 53 void defineImage( int nx=-1, … … 74 76 75 77 void getData( String &infile, 76 Matrix<Float> &spectra,78 Cube<Float> &spectra, 77 79 Matrix<Double> &direction, 78 Matrix<uChar> &flagtra,79 Vector<uInt> &rflag ) ;80 Cube<uChar> &flagtra, 81 Matrix<uInt> &rflag ) ; 80 82 81 void toInt( Matrix<uChar> &u, Matrix<Int> &v ) ;82 void toInt( Vector<uInt> &u, Vector<Int> &v ) ;83 void toInt( Array<uChar> *u, Array<Int> *v ) ; 84 void toInt( Array<uInt> *u, Array<Int> *v ) ; 83 85 84 86 void toPixel( Matrix<Double> &world, Matrix<Double> &pixel ) ; … … 101 103 Int convSampling_ ; 102 104 Array<Float> data_ ; 105 Vector<uInt> pollist_ ; 103 106 104 107 Table tab_ ; -
trunk/src/python_STGrid.cpp
r2356 r2360 3 3 //#--------------------------------------------------------------------------- 4 4 #include <string> 5 #include <vector> 5 6 6 7 #include <boost/python.hpp> … … 24 25 // .def("_setoption", &STGridWrapper::setOption) 25 26 // .def("_grid", &STGridWrapper::grid) 27 .def("_setpollist", &STGrid::setPolList) 26 28 .def("_defineimage", &STGrid::defineImage) 27 29 .def("_setoption", &STGrid::setOption)
Note:
See TracChangeset
for help on using the changeset viewer.