Changeset 2360


Ignore:
Timestamp:
12/06/11 13:13:24 (12 years ago)
Author:
Takeshi Nakazato
Message:

New Development: No

JIRA Issue: Yes CAS-2816

Ready for Test: No

Interface Changes: Yes

What Interface Changed: Added interface to select polarization

Test Programs: List test programs

Put in Release Notes: Yes/No?

Module(s): Module Names change impacts.

Description: Describe your changes here...

Polarization selection support.
Default is all polarization components.


Location:
trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/python/asapgrid.py

    r2358 r2360  
    1313        self.gridder._setin( infile )
    1414
     15    def setPolList( self, pollist ):
     16        self.gridder._setpollist( pollist )
     17
    1518    def defineImage( self, nx=-1, ny=-1, cellx='', celly='', center='' ):
    1619        self.gridder._defineimage( nx, ny, cellx, celly, center )
     
    2528        self.outfile = self.gridder._save( outfile )
    2629
    27     def plot( self, plotchan=-1 ):
     30    def plot( self, plotchan=-1, plotpol=-1 ):
    2831        plotter = _SDGridPlotter( self.infile, self.outfile )
    29         plotter.plot( chan=plotchan )
     32        plotter.plot( chan=plotchan, pol=plotpol )
    3033       
    3134class _SDGridPlotter:
     
    4144        self.ny = -1
    4245        self.nchan = 0
     46        self.npol = 0
     47        self.pollist = []
    4348        self.cellx = 0.0
    4449        self.celly = 0.0
     
    5762        del s
    5863
    59         idx = spectra.nonzero()[1]
    60         #self.nonzero = self.pointing.take( idx, axis=1 )
    61        
    6264        s = scantable( self.outfile, average=False )
    63         self.grid = numpy.array( s.get_directionval() ).transpose()
    64         dirstring = s.get_direction()
    6565        nrow = s.nrow()
    66         spectra = []
     66        pols = numpy.ones( nrow, dtype=int )
    6767        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
    7083
    7184        idx = 0
     
    7386        while ( dirstring[idx].split()[-1] == d0 ): 
    7487            idx += 1
    75 
    7688        self.ny = idx
    77         self.nx = nrow / idx
     89        self.nx = nrow / (self.npol * idx )
     90        #print 'nx,ny=',self.nx,self.ny
    7891       
    7992        self.cellx = abs( self.grid[0][0] - self.grid[0][1] )
    8093        self.celly = abs( self.grid[1][0] - self.grid[1][self.ny] )
     94        #print 'cellx,celly=',self.cellx,self.celly
    8195
    82         self.data = spectra.reshape( (self.nchan,self.nx,self.ny) )
     96        self.data = spectra.reshape( (self.npol,self.nchan,self.nx,self.ny) )
    8397
    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)
    85107        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'
    88110        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)
    91114        pl.figure(10)
    92115        pl.clf()
    93116        pl.plot(self.grid[0],self.grid[1],'.',color='blue')
    94117        pl.plot(self.pointing[0],self.pointing[1],'.',color='red')
    95         #pl.plot(self.nonzero[0],self.nonzero[1],'o',color='green')
    96118        extent=[self.grid[0].min()-0.5*self.cellx,
    97119                self.grid[0].max()+0.5*self.cellx,
  • trunk/src/STGrid.cpp

    r2359 r2360  
    1313
    1414#include <tables/Tables/Table.h>
     15#include <tables/Tables/TableRecord.h>
     16#include <tables/Tables/ExprNode.h>
    1517#include <tables/Tables/ScalarColumn.h>
    1618#include <tables/Tables/ArrayColumn.h>
     
    6062    tab_ = Table( infile_ ) ;
    6163  }
     64}
     65
     66void STGrid::setPolList( vector<unsigned int> pols )
     67{
     68  //pollist_ = Vector<uInt>( pols ) ;
     69  pollist_.assign( Vector<uInt>( pols ) ) ;
     70  cout << "pollist_ = " << pollist_ << endl ;
    6271}
    6372
     
    126135
    127136  // retrieve data
    128   Matrix<Float> spectra ;
     137  Cube<Float> spectra ;
    129138  Matrix<Double> direction ;
    130   Matrix<uChar> flagtra ;
    131   Vector<uInt> rflag ;
     139  Cube<uChar> flagtra ;
     140  Matrix<uInt> rflag ;
    132141  getData( infile_, spectra, direction, flagtra, rflag ) ;
    133142  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 ;
    137148
    138149  // flagtra: uChar -> Int
    139150  // 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 ) ;
    144157 
    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 ;
    147160 
    148161  // grid parameter
     
    177190  Matrix<Float> weight( IPosition( 2, nchan, nrow ), 1.0 ) ;
    178191
    179   // call ggridsd
     192//   // call ggridsd
    180193  Bool deletePos, deleteData, deleteWgt, deleteFlag, deleteFlagR, deleteConv, deleteDataG, deleteWgtG ;
    181194  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 ) ;
    183197  setReal( dataC, spectra ) ;
    184198  const Complex *data_p = dataC.getStorage( deleteData ) ;
     
    187201  const Int *rflag_p = rflagI.getStorage( deleteFlagR ) ;
    188202  Float *conv_p = convFunc.getStorage( deleteConv ) ;
    189   Int npol = 1 ;
     203  //Int npol = 1 ;
    190204  // Extend grid plane with convSupport_
    191205  //IPosition gshape( 4, nx_, ny_, npol, nchan ) ;
     
    235249           gdata_p,
    236250           wdata_p,
    237            //&nx_,
    238            //&ny_,
    239251           &gnx,
    240252           &gny,
     
    381393
    382394void STGrid::getData( String &infile,
    383                       Matrix<Float> &spectra,
     395                      Cube<Float> &spectra,
    384396                      Matrix<Double> &direction,
    385                       Matrix<uChar> &flagtra,
    386                       Vector<uInt> &rflag )
     397                      Cube<uChar> &flagtra,
     398                      Matrix<uInt> &rflag )
    387399{
    388400  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
     460void STGrid::toInt( Array<uChar> *u, Array<Int> *v )
     461{
     462  uInt len = u->nelements() ;
    402463  Int *int_p = new Int[len] ;
    403464  Bool deleteIt ;
    404   const uChar *data_p = u.getStorage( deleteIt ) ;
     465  const uChar *data_p = u->getStorage( deleteIt ) ;
    405466  Int *i_p = int_p ;
    406467  const uChar *u_p = data_p ;
     
    410471    u_p++ ;
    411472  }
    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
     477void STGrid::toInt( Array<uInt> *u, Array<Int> *v )
     478{
     479  uInt len = u->nelements() ;
    419480  Int *int_p = new Int[len] ;
    420481  Bool deleteIt ;
    421   const uInt *data_p = u.getStorage( deleteIt ) ;
     482  const uInt *data_p = u->getStorage( deleteIt ) ;
    422483  Int *i_p = int_p ;
    423484  const uInt *u_p = data_p ;
     
    427488    u_p++ ;
    428489  }
    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 ) ;
    431492}
    432493
     
    532593string STGrid::saveData( string outfile )
    533594{
    534   Int polno = 0 ;
     595  //Int polno = 0 ;
    535596  string outfile_ ;
    536597  if ( outfile.size() == 0 ) {
     
    550611  CountedPtr<Scantable> out( new Scantable( *ref, True ) ) ;
    551612  Table tab = out->table() ;
    552   Int nrow = nx_ * ny_ ;
    553   tab.addRow( nrow ) ;
    554613  IPosition dshape = data_.shape() ;
    555614  Int npol = dshape[2] ;
    556615  Int nchan = dshape[3] ;
     616  Int nrow = nx_ * ny_ * npol ;
     617  tab.rwKeywordSet().define( "nPol", npol ) ;
     618  tab.addRow( nrow ) ;
    557619  Vector<Double> cpix( 2 ) ;
    558620  cpix(0) = Double( nx_ - 1 ) * 0.5 ;
     
    574636        spectraCol.put( irow, sp ) ;
    575637        directionCol.put( irow, dir ) ;
    576         polnoCol.put( irow, polno ) ;
     638        polnoCol.put( irow, pollist_[ipol] ) ;
    577639        irow++ ;
    578640      }
  • trunk/src/STGrid.h

    r2356 r2360  
    1616#include <fstream>
    1717#include <string>
     18#include <vector>
    1819
    1920#include <casa/BasicSL/String.h>
    2021#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>
    2324// #include <casa/Arrays/ArrayMath.h>
    2425// #include <casa/Inputs/Input.h>
     
    4748
    4849  void setFileIn( const string infile ) ;
    49 //   void setFileOut( const string outfile ) ;
     50
     51  void setPolList( vector<unsigned int> pols ) ;
    5052
    5153  void defineImage( int nx=-1,
     
    7476 
    7577  void getData( String &infile,
    76                 Matrix<Float> &spectra,
     78                Cube<Float> &spectra,
    7779                Matrix<Double> &direction,
    78                 Matrix<uChar> &flagtra,
    79                 Vector<uInt> &rflag ) ;
     80                Cube<uChar> &flagtra,
     81                Matrix<uInt> &rflag ) ;
    8082
    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 ) ;
    8385
    8486  void toPixel( Matrix<Double> &world, Matrix<Double> &pixel ) ;
     
    101103  Int convSampling_ ;
    102104  Array<Float> data_ ;
     105  Vector<uInt> pollist_ ;
    103106
    104107  Table tab_ ;
  • trunk/src/python_STGrid.cpp

    r2356 r2360  
    33//#---------------------------------------------------------------------------
    44#include <string>
     5#include <vector>
    56
    67#include <boost/python.hpp>
     
    2425//     .def("_setoption", &STGridWrapper::setOption)
    2526//     .def("_grid", &STGridWrapper::grid)
     27    .def("_setpollist", &STGrid::setPolList)
    2628    .def("_defineimage", &STGrid::defineImage)
    2729    .def("_setoption", &STGrid::setOption)
Note: See TracChangeset for help on using the changeset viewer.