Changes in trunk/src/STGrid.cpp [2461:2371]

File:

• trunk/src/STGrid.cpp (modified) (29 diffs)

trunk/src/STGrid.cpp

@@ -1 @@
-//
-// C++ Implementation: STGrid
-//
-// Description:
-//
-//
-// Author: Takeshi Nakazato <takeshi.nakazato@nao.ac.jp>, (C) 2011
-//
-// Copyright: See COPYING file that comes with this distribution
-//
-//
+#include <iostream>
+#include <fstream>
+
 #include <casa/BasicSL/String.h>
 #include <casa/Arrays/Vector.h>
+#include <casa/Arrays/Matrix.h>
+#include <casa/Arrays/Cube.h>
 #include <casa/Arrays/ArrayMath.h>
+#include <casa/Arrays/ArrayPartMath.h>
 #include <casa/Quanta/Quantum.h>
 #include <casa/Quanta/QuantumHolder.h>
@@ -20 +15 @@
 #include <tables/Tables/Table.h>
 #include <tables/Tables/TableRecord.h>
-#include <tables/Tables/TableRow.h>
 #include <tables/Tables/ExprNode.h>
 #include <tables/Tables/ScalarColumn.h>
 #include <tables/Tables/ArrayColumn.h>
-#include <tables/Tables/TableCopy.h>
 
 #include <measures/Measures/MDirection.h>
 
-#include "MathUtils.h"
-#include <atnf/PKSIO/SrcType.h>
+#include <MathUtils.h>
 
 #include "STGrid.h"
 
 using namespace std ;
-using namespace concurrent ;
 using namespace casa ;
 using namespace asap ;
@@ -40 +31 @@
 namespace asap {
 
-// for performance check
-double eToInt = 0.0 ;
-double eGetWeight = 0.0 ;
-
 // constructor
 STGrid::STGrid()
-  : vshape_( 1 ), wshape_( 2 ), dshape_( 2 )
 {
   init() ;
@@ -52 +38 @@
 
 STGrid::STGrid( const string infile )
-  : vshape_( 1 ), wshape_( 2 ), dshape_( 2 )
 {
   init() ;
 
   setFileIn( infile ) ;
-}
-
-STGrid::STGrid( const vector<string> infile )
-{
-  init() ;
-
-  setFileList( infile ) ;
 }
 
@@ -82 +60 @@
   userSupport_ = -1 ;
   convSampling_ = 100 ;
-  nprocessed_ = 0 ;
-  nchunk_ = 0 ;
-
-  // initialize user input 
-  nxUI_ = -1 ;
-  nyUI_ = -1 ;
-  cellxUI_ = "" ;
-  cellyUI_ = "" ;
-  centerUI_ = "" ;
-  doclip_ = False ;
 }
 
 void STGrid::setFileIn( const string infile )
 {
-  nfile_ = 1 ;
   String name( infile ) ;
-  infileList_.resize( nfile_ ) ;
-  infileList_[0] = String(infile) ;
-}
-
-void STGrid::setFileList( const vector<string> infile )
-{
-  nfile_ = infile.size() ;
-  infileList_.resize( nfile_ ) ;
-  for ( uInt i = 0 ; i < nfile_ ; i++ ) {
-    infileList_[i] = infile[i] ;
+  if ( infile_.compare( name ) != 0 ) {
+    infile_ = String( infile ) ;
+    tab_ = Table( infile_ ) ;
   }
 }
@@ -114 +74 @@
 {
   pollist_.assign( Vector<uInt>( pols ) ) ;
+  cout << "pollist_ = " << pollist_ << endl ;
 }
 
@@ -119 +80 @@
 {
   scanlist_.assign( Vector<uInt>( scans ) ) ;
+  cout << "scanlist_ = " << scanlist_ << endl ;
 }
 
@@ -125 +87 @@
   wtype_ = String( wType ) ;
   wtype_.upcase() ;
+  cout << "wtype_ = " << wtype_ << endl ; 
 }
 
@@ -133 +96 @@
                           string scenter ) 
 {
-  nxUI_ = (Int)nx ;
-  nyUI_ = (Int)ny ;
-  cellxUI_ = String( scellx ) ;
-  cellyUI_ = String( scelly ) ;
-  centerUI_ = String( scenter ) ;
+  ROArrayColumn<Double> dirCol( tab_, "DIRECTION" ) ;
+  Matrix<Double> direction = dirCol.getColumn() ;
+  Double rmax, rmin, dmax, dmin ;
+  minMax( rmin, rmax, direction.row( 0 ) ) ;
+  minMax( dmin, dmax, direction.row( 1 ) ) ;
+
+  Int npx = (Int)nx ;
+  Int npy = (Int)ny ;
+  String cellx( scellx ) ;
+  String celly( scelly ) ;
+  String center( scenter ) ;
+  setupGrid( npx, npy, 
+             cellx, celly, 
+             rmin, rmax, 
+             dmin, dmax, 
+             center ) ;
 }
   
@@ -176 +150 @@
                 Double*);
 }
-void STGrid::call_ggridsd( Array<Double> &xypos,
-                           Array<Complex> &spectra,
-                           Int &nvispol,
-                           Int &nvischan,
-                           Array<Int> &flagtra,
-                           Array<Int> &flagrow,
-                           Array<Float> &weight,
-                           Int &nrow,
-                           Int &irow,
-                           Array<Complex> &gdata,
-                           Array<Float> &gwgt,
-                           Int &nx,
-                           Int &ny,
-                           Int &npol,
-                           Int &nchan,
-                           Int &support,
-                           Int &sampling,
-                           Vector<Float> &convFunc,
-                           Int *chanMap,
-                           Int *polMap ) 
-{
-  // parameters for gridding
-  Int idopsf = 0 ;
-  Int len = npol*nchan ;
-  Double *sumw_p = new Double[len] ;
-  {
-    Double *work_p = sumw_p ;
-    for ( Int i = 0 ; i < len ; i++ ) {
-      *work_p = 0.0 ;
-      work_p++ ;
-    }
-  }
-
-  // prepare pointer
-  Bool deletePos, deleteData, deleteWgt, deleteFlag, deleteFlagR, deleteConv, deleteDataG, deleteWgtG ;
-  Double *xy_p = xypos.getStorage( deletePos ) ;
-  const Complex *values_p = spectra.getStorage( deleteData ) ;
-  const Int *flag_p = flagtra.getStorage( deleteFlag ) ;
-  const Int *rflag_p = flagrow.getStorage( deleteFlagR ) ;
-  const Float *wgt_p = weight.getStorage( deleteWgt ) ;
-  Complex *grid_p = gdata.getStorage( deleteDataG ) ;
-  Float *wgrid_p = gwgt.getStorage( deleteWgtG ) ;
-  Float *conv_p = convFunc.getStorage( deleteConv ) ;
-
-  // pass copy of irow to ggridsd since it will be modified in theroutine
-  Int irowCopy = irow ;
-      
-  // call ggridsd
-  ggridsd( xy_p,
-           values_p,
-           &nvispol,
-           &nvischan,
-           &idopsf,
-           flag_p,
-           rflag_p,
-           wgt_p,
-           &nrow,
-           &irowCopy,
-           grid_p,
-           wgrid_p,
-           &nx,
-           &ny,
-           &npol,
-           &nchan,
-           &support,
-           &sampling,
-           conv_p,
-           chanMap,
-           polMap,
-           sumw_p ) ;
-
-  // finalization
-  xypos.putStorage( xy_p, deletePos ) ;
-  spectra.freeStorage( values_p, deleteData ) ;
-  flagtra.freeStorage( flag_p, deleteFlag ) ;
-  flagrow.freeStorage( rflag_p, deleteFlagR ) ;
-  weight.freeStorage( wgt_p, deleteWgt ) ;
-  gdata.putStorage( grid_p, deleteDataG ) ;
-  gwgt.putStorage( wgrid_p, deleteWgtG ) ;
-  convFunc.putStorage( conv_p, deleteConv ) ;
-  delete sumw_p ;
-}
-
-#define NEED_UNDERSCORES
-#if defined(NEED_UNDERSCORES)
-#define ggridsd2 ggridsd2_
-#endif
-extern "C" { 
-   void ggridsd2(Double*,
-                 const Complex*,
-                 Int*,
-                 Int*,
-                 Int*,
-                 const Int*,
-                 const Int*,
-                 const Float*,
-                 Int*,
-                 Int*,
-                 Complex*,
-                 Float*,
-                 Int*,
-                 Complex*,
-                 Float*,
-                 Float*,
-                 Complex*,
-                 Float*,
-                 Float*,
-                 Int*,
-                 Int*,
-                 Int *,
-                 Int *,
-                 Int*,
-                 Int*,
-                 Float*,
-                 Int*,
-                 Int*,
-                 Double*);
-}
-void STGrid::call_ggridsd2( Array<Double> &xypos,
-                            Array<Complex> &spectra,
-                            Int &nvispol,
-                            Int &nvischan,
-                            Array<Int> &flagtra,
-                            Array<Int> &flagrow,
-                            Array<Float> &weight,
-                            Int &nrow,
-                            Int &irow,
-                            Array<Complex> &gdata,
-                            Array<Float> &gwgt,
-                            Array<Int> &npoints,
-                            Array<Complex> &clipmin,
-                            Array<Float> &clipwmin,
-                            Array<Float> &clipcmin,
-                            Array<Complex> &clipmax,
-                            Array<Float> &clipwmax,
-                            Array<Float> &clipcmax,
-                            Int &nx,
-                            Int &ny,
-                            Int &npol,
-                            Int &nchan,
-                            Int &support,
-                            Int &sampling,
-                            Vector<Float> &convFunc,
-                            Int *chanMap,
-                            Int *polMap ) 
-{
-  // parameters for gridding
-  Int idopsf = 0 ;
-  Int len = npol*nchan ;
-  Double *sumw_p = new Double[len] ;
-  {
-    Double *work_p = sumw_p ;
-    for ( Int i = 0 ; i < len ; i++ ) {
-      *work_p = 0.0 ;
-      work_p++ ;
-    }
-  }
-
-  // prepare pointer
-  Bool deletePos, deleteData, deleteWgt, deleteFlag, deleteFlagR, deleteConv, deleteDataG, deleteWgtG, deleteNpts, deleteCMin, deleteCWMin, deleteCCMin, deleteCMax, deleteCWMax, deleteCCMax ;
-  Double *xy_p = xypos.getStorage( deletePos ) ;
-  const Complex *values_p = spectra.getStorage( deleteData ) ;
-  const Int *flag_p = flagtra.getStorage( deleteFlag ) ;
-  const Int *rflag_p = flagrow.getStorage( deleteFlagR ) ;
-  const Float *wgt_p = weight.getStorage( deleteWgt ) ;
-  Complex *grid_p = gdata.getStorage( deleteDataG ) ;
-  Float *wgrid_p = gwgt.getStorage( deleteWgtG ) ;
-  Float *conv_p = convFunc.getStorage( deleteConv ) ;
-  Int *npts_p = npoints.getStorage( deleteNpts ) ;
-  Complex *cmin_p = clipmin.getStorage( deleteCMin ) ;
-  Float *cwmin_p = clipwmin.getStorage( deleteCWMin ) ;
-  Float *ccmin_p = clipcmin.getStorage( deleteCCMin ) ;
-  Complex *cmax_p = clipmax.getStorage( deleteCMax ) ;
-  Float *cwmax_p = clipwmax.getStorage( deleteCWMax ) ;
-  Float *ccmax_p = clipcmax.getStorage( deleteCCMax ) ;
-
-  // pass copy of irow to ggridsd since it will be modified in theroutine
-  Int irowCopy = irow ;
-      
-  // call ggridsd
-  ggridsd2( xy_p,
-            values_p,
-            &nvispol,
-            &nvischan,
-            &idopsf,
-            flag_p,
-            rflag_p,
-            wgt_p,
-            &nrow,
-            &irowCopy,
-            grid_p,
-            wgrid_p,
-            npts_p,
-            cmin_p,
-            cwmin_p,
-            ccmin_p,
-            cmax_p,
-            cwmax_p,
-            ccmax_p,
-            &nx,
-            &ny,
-            &npol,
-            &nchan,
-            &support,
-            &sampling,
-            conv_p,
-            chanMap,
-            polMap,
-            sumw_p ) ;
-
-  // finalization
-  xypos.putStorage( xy_p, deletePos ) ;
-  spectra.freeStorage( values_p, deleteData ) ;
-  flagtra.freeStorage( flag_p, deleteFlag ) ;
-  flagrow.freeStorage( rflag_p, deleteFlagR ) ;
-  weight.freeStorage( wgt_p, deleteWgt ) ;
-  gdata.putStorage( grid_p, deleteDataG ) ;
-  gwgt.putStorage( wgrid_p, deleteWgtG ) ;
-  convFunc.putStorage( conv_p, deleteConv ) ;
-  clipmin.putStorage( cmin_p, deleteCMin ) ;
-  clipwmin.putStorage( cwmin_p, deleteCWMin ) ;
-  clipcmin.putStorage( ccmin_p, deleteCCMin ) ;
-  clipmax.putStorage( cmax_p, deleteCMax ) ;
-  clipwmax.putStorage( cwmax_p, deleteCWMax ) ;
-  clipcmax.putStorage( ccmax_p, deleteCCMax ) ;
-  delete sumw_p ;
-}
-
-void STGrid::grid()
+void STGrid::grid() 
 {
   LogIO os( LogOrigin("STGrid", "grid", WHERE) ) ;
-  double t0,t1 ;
-
-  // data selection
+
+  // retrieve data
+  Cube<Float> spectra ;
+  Matrix<Double> direction ;
+  Cube<uChar> flagtra ;
+  Matrix<uInt> rflag ;
+  Matrix<Float> weight ;
+  double t0, t1 ;
   t0 = mathutil::gettimeofday_sec() ;
-  selectData() ;
+  getData( spectra, direction, flagtra, rflag, weight ) ;
   t1 = mathutil::gettimeofday_sec() ;
-  os << LogIO::DEBUGGING << "selectData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ;
-
-  setupGrid() ;
-  setupArray() ;
-
-  if ( wtype_.compare("UNIFORM") != 0 &&
-       wtype_.compare("TINT") != 0 && 
-       wtype_.compare("TSYS") != 0 &&
-       wtype_.compare("TINTSYS") != 0 ) {
-    LogIO os( LogOrigin("STGrid", "grid", WHERE) ) ;
-    os << LogIO::WARN << "Unsupported weight type '" << wtype_ << "', apply UNIFORM weight" << LogIO::POST ;
-    wtype_ = "UNIFORM" ;
-  }
-
+  os << "getData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
+  IPosition sshape = spectra.shape() ;
+  //os << "spectra.shape()=" << spectra.shape() << LogIO::POST ;
+  //os << "max(spectra) = " << max(spectra) << LogIO::POST ;
+  //os << "weight = " << weight << LogIO::POST ;
+
+  // flagtra: uChar -> Int
+  // rflag: uInt -> Int
+  Cube<Int> flagI ;
+  Matrix<Int> rflagI ;
+  t0 = mathutil::gettimeofday_sec() ;
+  toInt( &flagtra, &flagI ) ;
+  toInt( &rflag, &rflagI ) ;
+  t1 = mathutil::gettimeofday_sec() ;
+  os << "toInt: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
+  
   // grid parameter
   os << LogIO::DEBUGGING ;
-  os << "----------" << endl ;
-  os << "Data selection summary" << endl ;
-  os << "   ifno = " << ifno_ << endl ;
-  os << "   pollist = " << pollist_ << endl ;
-  os << "   scanlist = " << scanlist_ << endl ;
   os << "----------" << endl ;
   os << "Grid parameter summary" << endl ;
@@ -439 +187 @@
   os << "   (cellx,celly) = (" << cellx_ << "," << celly_ << ")" << endl ;
   os << "   center = " << center_ << endl ;
-  os << "   weighting = " << wtype_ << endl ;
-  os << "   convfunc = " << convType_ << " with support " << convSupport_ << endl ; 
-  os << "   doclip = " << (doclip_?"True":"False") << endl ;
   os << "----------" << LogIO::POST ;
   os << LogIO::NORMAL ;
 
-  if ( doclip_ )
-    gridPerRowWithClipping() ;
-  else 
-    gridPerRow() ;
-}
-
-void STGrid::updateChunkShape()
-{
-  // TODO: nchunk_ must be determined from nchan_, npol_, and (nx_,ny_) 
-  //       by considering data size to be allocated for ggridsd input/output
-  nchunk_ = 400 ;
-  nchunk_ = min( nchunk_, nrow_ ) ;
-  vshape_ = IPosition( 1, nchunk_ ) ;
-  wshape_ = IPosition( 2, nchan_, nchunk_ ) ;
-  dshape_ = IPosition( 2, 2, nchunk_ ) ;
-}
-
-struct STGChunk {
-  Int nrow ;
-  Array<Complex> spectra;
-  Array<Int> flagtra;
-  Array<Int> rflag;
-  Array<Float> weight;
-  Array<Double> direction;
-  STGChunk(IPosition const &wshape, IPosition const &vshape,
-           IPosition const &dshape)
-    : spectra(wshape), flagtra(wshape), rflag(vshape), weight(wshape),
-      direction(dshape)
-  { }
-};
-
-struct STCommonData {
-  Int gnx;
-  Int gny;
-  Int *chanMap;
+  // convolution kernel
   Vector<Float> convFunc ;
-  Array<Complex> gdataArrC;
-  Array<Float> gwgtArr;
-  STCommonData(IPosition const &gshape, Array<Float> const &data)
-    : gdataArrC(gshape, 0.0), gwgtArr(data) {}
-};
-
-struct STCommonDataWithClipping {
-  Int gnx;
-  Int gny;
-  Int *chanMap;
-  Vector<Float> convFunc ;
-  Array<Complex> gdataArrC;
-  Array<Float> gwgtArr;
-  Array<Int> npoints ;
-  Array<Complex> clipMin ;
-  Array<Float> clipWMin ;
-  Array<Float> clipCMin ;
-  Array<Complex> clipMax ;
-  Array<Float> clipWMax ;
-  Array<Float> clipCMax ;  
-  STCommonDataWithClipping(IPosition const &gshape, 
-                           IPosition const &pshape, 
-                           Array<Float> const &data)
-    : gdataArrC(gshape, 0.0), 
-      gwgtArr(data), 
-      npoints(pshape, 0),
-      clipMin(gshape, Complex(FLT_MAX,0.0)),
-      clipWMin(gshape, 0.0),
-      clipCMin(gshape, 0.0),
-      clipMax(gshape, Complex(-FLT_MAX,0.0)),
-      clipWMax(gshape, 0.0),
-      clipCMax(gshape, 0.0)
-  {}
-};
-
-#define DO_AHEAD 3
-
-struct STContext {
-  STCommonData &common;
-  FIFO<STGChunk *, DO_AHEAD> queue;
-  STGrid *const self;
-  const Int pol;
-  STContext(STGrid *obj, STCommonData &common, Int pol)
-    : self(obj), common(common), pol(pol) {}
-};
-
-struct STContextWithClipping {
-  STCommonDataWithClipping &common;
-  FIFO<STGChunk *, DO_AHEAD> queue;
-  STGrid *const self;
-  const Int pol;
-  STContextWithClipping(STGrid *obj, STCommonDataWithClipping &common, Int pol)
-    : self(obj), common(common), pol(pol) {}
-};
-
-
-bool STGrid::produceChunk(void *ctx) throw(PCException)
-{
-  STContext &context = *(STContext *)ctx;
-  if ( context.self->nprocessed_ >= context.self->nrow_ ) {
-    return false;
-  }
-  STGChunk *chunk = new STGChunk(context.self->wshape_,
-                                 context.self->vshape_,
-                                 context.self->dshape_);
-
-  double t0 = mathutil::gettimeofday_sec() ;
-  chunk->nrow = context.self->getDataChunk(
-        context.self->wshape_, context.self->vshape_, context.self->dshape_,
-        chunk->spectra, chunk->direction,
-        chunk->flagtra, chunk->rflag, chunk->weight);
-  double t1 = mathutil::gettimeofday_sec() ;
-  context.self->eGetData_ += t1-t0 ;
-
-  context.queue.lock();
-  context.queue.put(chunk);
-  context.queue.unlock();
-  return true;
-}
-
-void STGrid::consumeChunk(void *ctx) throw(PCException)
-{
-  STContext &context = *(STContext *)ctx;
-  STGChunk *chunk = NULL;
-  try {
-    context.queue.lock();
-    chunk = context.queue.get();
-    context.queue.unlock();
-  } catch (FullException &e) {
-    context.queue.unlock();
-    // TODO: log error
-    throw PCException();
-  }
-
-  double t0, t1 ;
+  t0 = mathutil::gettimeofday_sec() ;
+  setConvFunc( convFunc ) ;
+  t1 = mathutil::gettimeofday_sec() ;
+  os << "setConvFunc: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
+  //cout << "convSupport=" << convSupport_ << endl ;
+  //cout << "convFunc=" << convFunc << endl ;
+
   // world -> pixel
-  Array<Double> xypos( context.self->dshape_ ) ;
+  Matrix<Double> xypos( direction.shape(), 0.0 ) ;
   t0 = mathutil::gettimeofday_sec() ;
-  context.self->toPixel( chunk->direction, xypos ) ;
+  toPixel( direction, xypos ) ;  
   t1 = mathutil::gettimeofday_sec() ;
-  context.self->eToPixel_ += t1-t0 ;
-   
+  os << "toPixel: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
+  
   // call ggridsd
-  Int nvispol = 1 ;
+  Bool deletePos, deleteData, deleteWgt, deleteFlag, deleteFlagR, deleteConv, deleteDataG, deleteWgtG ;
+  Double *xypos_p = xypos.getStorage( deletePos ) ;
+  Cube<Complex> dataC( spectra.shape(), 0.0 ) ;
+  setReal( dataC, spectra ) ;
+  const Complex *data_p = dataC.getStorage( deleteData ) ;
+  const Float *wgt_p = weight.getStorage( deleteWgt ) ;
+  const Int *flag_p = flagI.getStorage( deleteFlag ) ;
+  const Int *rflag_p = rflagI.getStorage( deleteFlagR ) ;
+  Float *conv_p = convFunc.getStorage( deleteConv ) ;
+  // Extend grid plane with convSupport_
+  Int gnx = nx_ ;
+  Int gny = ny_ ;
+//   Int gnx = nx_+convSupport_*2 ;
+//   Int gny = ny_+convSupport_*2 ;
+  IPosition gshape( 4, gnx, gny, npol_, nchan_ ) ;
+  Array<Complex> gdataArrC( gshape, 0.0 ) ;
+  Array<Float> gwgtArr( gshape, 0.0 ) ;
+  Complex *gdata_p = gdataArrC.getStorage( deleteDataG ) ;
+  Float *wdata_p = gwgtArr.getStorage( deleteWgtG ) ;
+  Int idopsf = 0 ;
+  Int *chanMap = new Int[nchan_] ;
+  {
+    Int *work_p = chanMap ;
+    for ( Int i = 0 ; i < nchan_ ; i++ ) {
+      *work_p = i ;
+      work_p++ ;
+    }
+  }
+  Int *polMap = new Int[npol_] ;
+  {
+    Int *work_p = polMap ;
+    for ( Int i = 0 ; i < npol_ ; i++ ) {
+      *work_p = i ;
+      work_p++ ;
+    }
+  }
+  Double *sumw_p = new Double[npol_*nchan_] ;
+  {
+    Double *work_p = sumw_p ;
+    for ( Int i = 0 ; i < npol_*nchan_ ; i++ ) {
+      *work_p = 0.0 ;
+      work_p++ ;
+    }
+  }
+  t0 = mathutil::gettimeofday_sec() ;
   Int irow = -1 ;
-  t0 = mathutil::gettimeofday_sec() ;
-  context.self->call_ggridsd( xypos,
-                chunk->spectra,
-                nvispol,
-                context.self->nchan_,
-                chunk->flagtra,
-                chunk->rflag,
-                chunk->weight,
-                chunk->nrow,
-                irow,
-                context.common.gdataArrC,
-                context.common.gwgtArr,
-                context.common.gnx,
-                context.common.gny,
-                context.self->npol_,
-                context.self->nchan_,
-                context.self->convSupport_,
-                context.self->convSampling_,
-                context.common.convFunc,
-                context.common.chanMap,
-                (Int*)&context.pol ) ;
+  ggridsd( xypos_p,
+           data_p,
+           &npol_,
+           &nchan_,
+           &idopsf,
+           flag_p,
+           rflag_p,
+           wgt_p,
+           &nrow_,
+           &irow,
+           gdata_p,
+           wdata_p, 
+           &gnx,
+           &gny,
+           &npol_,
+           &nchan_,
+           &convSupport_,
+           &convSampling_,
+           conv_p,
+           chanMap,
+           polMap,
+           sumw_p ) ;
   t1 = mathutil::gettimeofday_sec() ;
-  context.self->eGGridSD_ += t1-t0 ;
-  
-  delete chunk;
-}
-
-void STGrid::gridPerRow()
-{
-  LogIO os( LogOrigin("STGrid", "gridPerRow", WHERE) ) ;
-  double t0, t1 ;
-
-
-  // grid data
-  // Extend grid plane with convSupport_
-  //   Int gnx = nx_+convSupport_*2 ;
-  //   Int gny = ny_+convSupport_*2 ;
-  Int gnx = nx_;
-  Int gny = ny_;
-
-  IPosition gshape( 4, gnx, gny, npol_, nchan_ ) ;
-  // 2011/12/20 TN
-  // data_ and gwgtArr share storage
-  data_.resize( gshape ) ;
-  data_ = 0.0 ;
-  STCommonData common = STCommonData(gshape, data_);
-  common.gnx = gnx ;
-  common.gny = gny ;
-
-  // parameters for gridding
-  Int *chanMap = new Int[nchan_] ;
-  for ( Int i = 0 ; i < nchan_ ; i++ ) {
-    chanMap[i] = i ;
-  }
-  common.chanMap = chanMap;
-
-  // convolution kernel
-  t0 = mathutil::gettimeofday_sec() ;
-  setConvFunc( common.convFunc ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  os << LogIO::DEBUGGING << "setConvFunc: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
-
-  // for performance check
-  eGetData_ = 0.0 ;
-  eToPixel_ = 0.0 ;
-  eGGridSD_ = 0.0 ;
-  double eInitPol = 0.0 ;
-
-  for ( uInt ifile = 0 ; ifile < nfile_ ; ifile++ ) {
-    initTable( ifile ) ;
-
-    os << "start table " << ifile << ": " << infileList_[ifile] << LogIO::POST ;   
-    Broker broker = Broker(produceChunk, consumeChunk);
-    for ( Int ipol = 0 ; ipol < npol_ ; ipol++ ) {
-      t0 = mathutil::gettimeofday_sec() ;
-      initPol( ipol ) ; // set ptab_ and attach()
-      t1 = mathutil::gettimeofday_sec() ;
-      eInitPol += t1-t0 ;
-      
-      STContext context(this, common, ipol);
-      
-      os << "start pol " << ipol << LogIO::POST ;
-      
-      nprocessed_ = 0 ;
-#if 1
-      broker.runProducerAsMasterThread(&context, DO_AHEAD);
-#else
-      for (;;) {
-        bool produced = produceChunk(&context);
-        if (! produced) {
-          break;
-        }
-        consumeChunk(&context);
-      }
-#endif
-
-      os << "end pol " << ipol << LogIO::POST ;
-
-    }
-    os << "end table " << ifile << LogIO::POST ;   
-  }
-  os << LogIO::DEBUGGING << "initPol: elapsed time is " << eInitPol << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "getData: elapsed time is " << eGetData_-eToInt-eGetWeight << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "toPixel: elapsed time is " << eToPixel_ << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "ggridsd: elapsed time is " << eGGridSD_ << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "toInt: elapsed time is " << eToInt << " sec." << LogIO::POST ;
-  os << LogIO::DEBUGGING << "getWeight: elapsed time is " << eGetWeight << " sec." << LogIO::POST ;
-  
+  os << "ggridsd: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
+  xypos.putStorage( xypos_p, deletePos ) ;
+  dataC.freeStorage( data_p, deleteData ) ;
+  weight.freeStorage( wgt_p, deleteWgt ) ;
+  flagI.freeStorage( flag_p, deleteFlag ) ;
+  rflagI.freeStorage( rflag_p, deleteFlagR ) ;
+  convFunc.putStorage( conv_p, deleteConv ) ;
+  delete polMap ;
   delete chanMap ;
-
-  // set data
-  setData( common.gdataArrC, common.gwgtArr ) ;
-
-}
-
-void STGrid::consumeChunkWithClipping(void *ctx) throw(PCException)
-{
-  STContextWithClipping &context = *(STContextWithClipping *)ctx;
-  STGChunk *chunk = NULL;
-  try {
-    context.queue.lock();
-    chunk = context.queue.get();
-    context.queue.unlock();
-  } catch (FullException &e) {
-    context.queue.unlock();
-    // TODO: log error
-    throw PCException();
-  }
-
-  double t0, t1 ;
-  // world -> pixel
-  Array<Double> xypos( context.self->dshape_ ) ;
-  t0 = mathutil::gettimeofday_sec() ;
-  context.self->toPixel( chunk->direction, xypos ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  context.self->eToPixel_ += t1-t0 ;
-   
-  // call ggridsd
-  Int nvispol = 1 ;
-  Int irow = -1 ;
-  t0 = mathutil::gettimeofday_sec() ;
-  context.self->call_ggridsd2( xypos,
-                chunk->spectra,
-                nvispol,
-                context.self->nchan_,
-                chunk->flagtra,
-                chunk->rflag,
-                chunk->weight,
-                chunk->nrow,
-                irow,
-                context.common.gdataArrC,
-                context.common.gwgtArr,
-                context.common.npoints,
-                context.common.clipMin,
-                context.common.clipWMin,
-                context.common.clipCMin,
-                context.common.clipMax,
-                context.common.clipWMax,
-                context.common.clipCMax,
-                context.common.gnx,
-                context.common.gny,
-                context.self->npol_,
-                context.self->nchan_,
-                context.self->convSupport_,
-                context.self->convSampling_,
-                context.common.convFunc,
-                context.common.chanMap,
-                (Int*)&context.pol ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  context.self->eGGridSD_ += t1-t0 ;
-  
-  delete chunk;
-}
-
-void STGrid::gridPerRowWithClipping()
-{
-  LogIO os( LogOrigin("STGrid", "gridPerRowWithClipping", WHERE) ) ;
-  double t0, t1 ;
-
-
-  // grid data
-  // Extend grid plane with convSupport_
-  //   Int gnx = nx_+convSupport_*2 ;
-  //   Int gny = ny_+convSupport_*2 ;
-  Int gnx = nx_;
-  Int gny = ny_;
-
-  IPosition gshape( 4, gnx, gny, npol_, nchan_ ) ;
-  IPosition pshape( 3, gnx, gny, npol_ ) ;
-  // 2011/12/20 TN
-  // data_ and gwgtArr share storage
-  data_.resize( gshape ) ;
-  data_ = 0.0 ;
-  STCommonDataWithClipping common = STCommonDataWithClipping( gshape,
-                                                              pshape, 
-                                                              data_ ) ;
-  common.gnx = gnx ;
-  common.gny = gny ;
-
-  // parameters for gridding
-  Int *chanMap = new Int[nchan_] ;
-  for ( Int i = 0 ; i < nchan_ ; i++ ) {
-    chanMap[i] = i ;
-  }
-  common.chanMap = chanMap;
-
-  // convolution kernel
-  t0 = mathutil::gettimeofday_sec() ;
-  setConvFunc( common.convFunc ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  os << LogIO::DEBUGGING << "setConvFunc: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
-
-  // for performance check
-  eGetData_ = 0.0 ;
-  eToPixel_ = 0.0 ;
-  eGGridSD_ = 0.0 ;
-  double eInitPol = 0.0 ;
-
-  for ( uInt ifile = 0 ; ifile < nfile_ ; ifile++ ) {
-    initTable( ifile ) ;
-
-    os << "start table " << ifile << ": " << infileList_[ifile] << LogIO::POST ;   
-    Broker broker = Broker(produceChunk, consumeChunkWithClipping);
-    for ( Int ipol = 0 ; ipol < npol_ ; ipol++ ) {
-      t0 = mathutil::gettimeofday_sec() ;
-      initPol( ipol ) ; // set ptab_ and attach()
-      t1 = mathutil::gettimeofday_sec() ;
-      eInitPol += t1-t0 ;
-      
-      STContextWithClipping context(this, common, ipol);
-      
-      os << "start pol " << ipol << LogIO::POST ;
-      
-      nprocessed_ = 0 ;
-#if 1
-      broker.runProducerAsMasterThread(&context, DO_AHEAD);
-#else
-      for (;;) {
-        bool produced = produceChunk(&context);
-        if (! produced) {
-          break;
-        }
-        consumeChunkWithClipping(&context);
-      }
-#endif
-
-      os << "end pol " << ipol << LogIO::POST ;
-
-    }
-    os << "end table " << ifile << LogIO::POST ;   
-  }
-  os << LogIO::DEBUGGING << "initPol: elapsed time is " << eInitPol << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "getData: elapsed time is " << eGetData_-eToInt-eGetWeight << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "toPixel: elapsed time is " << eToPixel_ << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "ggridsd2: elapsed time is " << eGGridSD_ << " sec." << LogIO::POST ; 
-  os << LogIO::DEBUGGING << "toInt: elapsed time is " << eToInt << " sec." << LogIO::POST ;
-  os << LogIO::DEBUGGING << "getWeight: elapsed time is " << eGetWeight << " sec." << LogIO::POST ;
-  
-  delete chanMap ;
-
-  // clip min and max in each grid
-//   os << "BEFORE CLIPPING" << LogIO::POST ;
-//   os << "gdataArrC=" << common.gdataArrC << LogIO::POST ;
-//   os << "gwgtArr=" << common.gwgtArr << LogIO::POST ;
-  t0 = mathutil::gettimeofday_sec() ;
-  clipMinMax( common.gdataArrC, 
-              common.gwgtArr,
-              common.npoints,
-              common.clipMin,
-              common.clipWMin,
-              common.clipCMin,
-              common.clipMax,
-              common.clipWMax,
-              common.clipCMax ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  os << LogIO::DEBUGGING << "clipMinMax: elapsed time is " << t1-t0 << " sec." << LogIO::POST ;
-//   os << "AFTER CLIPPING" << LogIO::POST ;
-//   os << "gdataArrC=" << common.gdataArrC << LogIO::POST ;
-//   os << "gwgtArr=" << common.gwgtArr << LogIO::POST ;
-
-  // set data
-  setData( common.gdataArrC, common.gwgtArr ) ;
-
-}
-
-void STGrid::clipMinMax( Array<Complex> &grid,
-                         Array<Float> &weight,
-                         Array<Int> &npoints,
-                         Array<Complex> &clipmin,
-                         Array<Float> &clipwmin,
-                         Array<Float> &clipcmin,
-                         Array<Complex> &clipmax,
-                         Array<Float> &clipwmax,
-                         Array<Float> &clipcmax )
-{
-  //LogIO os( LogOrigin("STGrid","clipMinMax",WHERE) ) ;
-
-  // prepare pointers
-  Bool delG, delW, delNP, delCMin, delCWMin, delCCMin, delCMax, delCWMax, delCCMax ;
-  Complex *grid_p = grid.getStorage( delG ) ;
-  Float *wgt_p = weight.getStorage( delW ) ;
-  const Int *npts_p = npoints.getStorage( delNP ) ;
-  const Complex *cmin_p = clipmin.getStorage( delCMin ) ;
-  const Float *cwmin_p = clipwmin.getStorage( delCWMin ) ;
-  const Float *ccmin_p = clipcmin.getStorage( delCCMin ) ;
-  const Complex *cmax_p = clipmax.getStorage( delCMax ) ;
-  const Float *cwmax_p = clipwmax.getStorage( delCWMax ) ;
-  const Float *ccmax_p = clipcmax.getStorage( delCCMax ) ;
-
-  const IPosition &gshape = grid.shape() ;
-  long offset = gshape[0] * gshape[1] * gshape[2] ; // nx * ny * npol
-  Int nchan = gshape[3] ;
-  long origin = nchan * offset ;
-  for ( long i = 0 ; i < offset ; i++ ) {
-    if ( *npts_p > 2 ) {
-      for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
-        // clip minimum and maximum
-        *grid_p -= (*cmin_p)*(*cwmin_p)*(*ccmin_p)
-          + (*cmax_p)*(*cwmax_p)*(*ccmax_p) ;
-        *wgt_p -= (*cwmin_p)*(*ccmin_p)
-          + (*cwmax_p)*(*ccmax_p) ;
-        
-        grid_p += offset ;
-        wgt_p += offset ;
-        cmin_p += offset ;
-        cwmin_p += offset ;
-        ccmin_p += offset ;
-        cmax_p += offset ;
-        cwmax_p += offset ;
-        ccmax_p += offset ;
-      }
-      grid_p -= origin ;
-      wgt_p -= origin ;
-      cmin_p -= origin ;
-      cwmin_p -= origin ;
-      ccmin_p -= origin ;
-      cmax_p -= origin ;
-      cwmax_p -= origin ;
-      ccmax_p -= origin ;
-    }
-    grid_p++ ;
-    wgt_p++ ;
-    npts_p++ ;
-    cmin_p++ ;
-    cwmin_p++ ;
-    ccmin_p++ ;
-    cmax_p++ ;
-    cwmax_p++ ;
-    ccmax_p++ ;
-  }
-  grid_p -= offset ;
-  wgt_p -= offset ;
-  npts_p -= offset ;
-  cmin_p -= offset ;
-  cwmin_p -= offset ;
-  ccmin_p -= offset ;
-  cmax_p -= offset ;
-  cwmax_p -= offset ;
-  ccmax_p -= offset ;  
-
-  // finalization
-  grid.putStorage( grid_p, delG ) ;
-  weight.putStorage( wgt_p, delW ) ;
-  npoints.freeStorage( npts_p, delNP ) ;
-  clipmin.freeStorage( cmin_p, delCMin ) ;
-  clipwmin.freeStorage( cwmin_p, delCWMin ) ;
-  clipcmin.freeStorage( ccmin_p, delCCMin ) ;
-  clipmax.freeStorage( cmax_p, delCMax ) ;
-  clipwmax.freeStorage( cwmax_p, delCWMax ) ;
-  clipcmax.freeStorage( ccmax_p, delCCMax ) ;
-}
-
-void STGrid::initPol( Int ipol ) 
-{
-  LogIO os( LogOrigin("STGrid","initPol",WHERE) ) ;
-  if ( npolOrg_ == 1 ) {
-    os << "single polarization data." << LogIO::POST ;
-    ptab_ = tab_ ;
-  }
-  else 
-    ptab_ = tab_( tab_.col("POLNO") == pollist_[ipol] ) ;
-
-  attach( ptab_ ) ;
-}
-
-void STGrid::initTable( uInt idx ) 
-{
-  tab_ = tableList_[idx] ;
-  nrow_ = rows_[idx] ;
-  updateChunkShape() ;
-}
-
-void STGrid::setData( Array<Complex> &gdata,
+  gdataArrC.putStorage( gdata_p, deleteDataG ) ;
+  gwgtArr.putStorage( wdata_p, deleteWgtG ) ;
+  Array<Float> gdataArr = real( gdataArrC ) ;
+  setData( data_, gdataArr, gwgtArr ) ;
+  //Matrix<Double> sumWeight( IPosition( 2, npol_, nchan_ ), sumw_p, TAKE_OVER ) ;
+  delete sumw_p ;
+  //cout << "sumWeight = " << sumWeight << endl ;
+//   os << "gdataArr = " << gdataArr << LogIO::POST ;
+//   os << "gwgtArr = " << gwgtArr << LogIO::POST ;
+//   os << "data_ " << data_ << LogIO::POST ;
+}
+
+void STGrid::setData( Array<Float> &data,
+                      Array<Float> &gdata,
                       Array<Float> &gwgt )
 {
-  // 2011/12/20 TN
-  // gwgt and data_ share storage
   LogIO os( LogOrigin("STGrid","setData",WHERE) ) ;
   double t0, t1 ;
   t0 = mathutil::gettimeofday_sec() ;
-  uInt len = data_.nelements() ;
-  const Complex *w1_p ;
-  Float *w2_p ;
-  Bool b1, b2 ;
-  const Complex *gdata_p = gdata.getStorage( b1 ) ;
-  Float *gwgt_p = data_.getStorage( b2 ) ;
+  data.resize( gdata.shape() ) ;
+  uInt len = data.nelements() ;
+  Float *w0_p ;
+  const Float *w1_p, *w2_p ;
+  Bool b0, b1, b2 ;
+  Float *data_p = data.getStorage( b0 ) ;
+  const Float *gdata_p = gdata.getStorage( b1 ) ;
+  const Float *gwgt_p = gwgt.getStorage( b2 ) ;
+  w0_p = data_p ;
   w1_p = gdata_p ;
   w2_p = gwgt_p ;
   for ( uInt i = 0 ; i < len ; i++ ) {
-    if ( *w2_p > 0.0 ) *w2_p = (*w1_p).real() / *w2_p ;
+    *w0_p = (*w2_p > 0.0) ? (*w1_p / *w2_p) : 0.0 ;
+    w0_p++ ;
     w1_p++ ;
     w2_p++ ;
   }
+  data.putStorage( data_p, b0 ) ;
   gdata.freeStorage( gdata_p, b1 ) ;
-  data_.putStorage( gwgt_p, b2 ) ;
+  gwgt.freeStorage( gwgt_p, b2 ) ;
   t1 = mathutil::gettimeofday_sec() ;
-  os << LogIO::DEBUGGING << "setData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
-}
-
-void STGrid::setupGrid() 
-{
-  Double xmin,xmax,ymin,ymax ;
-  mapExtent( xmin, xmax, ymin, ymax ) ;
-  
-  setupGrid( nxUI_, nyUI_, cellxUI_, cellyUI_, 
-             xmin, xmax, ymin, ymax, centerUI_ ) ;
+  os << "setData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
 }
 
@@ -1024 +338 @@
                         String &center )
 {
-  LogIO os( LogOrigin("STGrid","setupGrid",WHERE) ) ;
   //cout << "nx=" << nx << ", ny=" << ny << endl ;
 
@@ -1058 +371 @@
     center_(0) = xcen.getValue( "rad" ) ;
     center_(1) = ycen.getValue( "rad" ) ;
-    double base = 0.5 * (xmin + xmax) ;
-    int maxrotate = 1 ;
-    int nelem = 2 * maxrotate + 1 ;
-    double *sep = new double[nelem] ;
-    for ( int i = 0 ; i < nelem ; i++ )
-      sep[i] = abs(base - center_[0] - (i-maxrotate) * C::_2pi) ;
-//     os << "sep[0]=" << sep[0] << endl  
-//        << "sep[1]=" << sep[1] << endl
-//        << "sep[2]=" << sep[2] << LogIO::POST ;
-    int idx = 0 ;
-    base = sep[0] ;
-    int nrotate = 0 ;
-    while ( idx < nelem ) {
-      if ( base > sep[idx] ) {
-        base = sep[idx] ;
-        nrotate = idx ;
-      }
-      idx++ ;
-    }
-    delete sep ;
-    nrotate -= maxrotate ;
-//     os << "nrotate = " << nrotate << LogIO::POST ;
-    center_[0] += nrotate * C::_2pi ;
-  }
-//   os << "xmin=" << xmin << LogIO::POST ;
-//   os << "center_=" << center_ << LogIO::POST ;
-
-
-  nx_ = nx ;
-  ny_ = ny ;
-  if ( nx < 0 && ny > 0 ) {
-    nx_ = ny ;
-    ny_ = ny ;
-  }
-  if ( ny < 0 && nx > 0 ) {
-    nx_ = nx ;
-    ny_ = nx ;
-  }
+  }
+
 
   //Double wx = xmax - xmin ;
@@ -1104 +381 @@
   wx *= 1.10 ;
   wy *= 1.10 ;
-
   Quantum<Double> qcellx ;
   Quantum<Double> qcelly ;
+  nx_ = nx ;
+  ny_ = ny ;
+  if ( nx < 0 && ny > 0 ) {
+    nx_ = ny ;
+    ny_ = ny ;
+  }
+  if ( ny < 0 && nx > 0 ) {
+    nx_ = nx ;
+    ny_ = nx ;
+  }
   //cout << "nx_ = " << nx_ << ",  ny_ = " << ny_ << endl ;
   if ( cellx.size() != 0 && celly.size() != 0 ) {
@@ -1113 +399 @@
   }
   else if ( celly.size() != 0 ) {
-    os << "Using celly to x-axis..." << LogIO::POST ;
+    cout << "Using celly to x-axis..." << endl ;
     readQuantity( qcelly, celly ) ;
     qcellx = qcelly ;
   }
   else if ( cellx.size() != 0 ) {
-    os << "Using cellx to y-axis..." << LogIO::POST ;
+    cout << "Using cellx to y-axis..." << endl ;
     readQuantity( qcellx, cellx ) ;
     qcelly = qcellx ;
@@ -1124 +410 @@
   else {
     if ( nx_ < 0 ) {
-      os << "No user preference in grid setting. Using default..." << LogIO::POST ;
+      cout << "No user preference in grid setting. Using default..." << endl ;
       readQuantity( qcellx, "1.0arcmin" ) ;
       qcelly = qcellx ;
     }
     else {
-      if ( wx == 0.0 ) {
-        os << "Using default spatial extent (10arcmin) in x" << LogIO::POST ;
-        wx = 0.00290888 ;
-      }
-      if ( wy == 0.0 ) {
-        os << "Using default spatial extent (10arcmin) in y" << LogIO::POST ;
-        wy = 0.00290888 ;
-      }
       qcellx = Quantum<Double>( wx/nx_, "rad" ) ;
       qcelly = Quantum<Double>( wy/ny_, "rad" ) ;
@@ -1142 +420 @@
   }
   cellx_ = qcellx.getValue( "rad" ) ;
-  // DEC correction
-  cellx_ /= cos( center_[1] ) ;
   celly_ = qcelly.getValue( "rad" ) ;
-  //os << "cellx_=" << cellx_ << ", celly_=" << celly_ << ", cos("<<center_(1)<<")=" << cos(center_(1)) << LogIO::POST ;
   if ( nx_ < 0 ) {
-    if ( wx == 0.0 ) {
-      os << "Using default spatial extent (10arcmin) in x" << LogIO::POST ;
-      wx = 0.00290888 ;
-    }
-    if ( wy == 0.0 ) {
-      os << "Using default spatial extent (10arcmin) in y" << LogIO::POST ;
-      wy = 0.00290888 ;
-    }
     nx_ = Int( ceil( wx/cellx_ ) ) ;
     ny_ = Int( ceil( wy/celly_ ) ) ;
@@ -1160 +427 @@
 }
 
-void STGrid::mapExtent( Double &xmin, Double &xmax, 
-                        Double &ymin, Double &ymax ) 
-{
-  //LogIO os( LogOrigin("STGrid","mapExtent",WHERE) ) ;
-  directionCol_.attach( tableList_[0], "DIRECTION" ) ;
-  Matrix<Double> direction = directionCol_.getColumn() ;
-  //os << "dirCol.nrow() = " << dirCol.nrow() << LogIO::POST ;
-  minMax( xmin, xmax, direction.row( 0 ) ) ;
-  minMax( ymin, ymax, direction.row( 1 ) ) ;
-  Double amin, amax, bmin, bmax ;
-  for ( uInt i = 1 ; i < nfile_ ; i++ ) {
-    directionCol_.attach( tableList_[i], "DIRECTION" ) ;
-    direction.assign( directionCol_.getColumn() ) ;
-    //os << "dirCol.nrow() = " << dirCol.nrow() << LogIO::POST ;
-    minMax( amin, amax, direction.row( 0 ) ) ;
-    minMax( bmin, bmax, direction.row( 1 ) ) ;
-    xmin = min( xmin, amin ) ;
-    xmax = max( xmax, amax ) ;
-    ymin = min( ymin, bmin ) ;
-    ymax = max( ymax, bmax ) ;
-  }
-  //os << "(xmin,xmax)=(" << xmin << "," << xmax << ")" << LogIO::POST ; 
-  //os << "(ymin,ymax)=(" << ymin << "," << ymax << ")" << LogIO::POST ; 
-}
-
-void STGrid::selectData()
-{
-  LogIO os( LogOrigin("STGrid","selectData",WHERE) ) ;    
+void STGrid::selectData( Table &tab )
+{
   Int ifno = ifno_ ;
-  tableList_.resize( nfile_ ) ;
-  if ( ifno_ == -1 ) {
-    Table taborg( infileList_[0] ) ;
+  Table taborg( infile_ ) ;
+  if ( ifno == -1 ) {
+    LogIO os( LogOrigin("STGrid","selectData",WHERE) ) ;
+//     os << LogIO::SEVERE
+//        << "Please set IFNO before actual gridding" 
+//        << LogIO::EXCEPTION ;
     ROScalarColumn<uInt> ifnoCol( taborg, "IFNO" ) ;
-    ifno_ = ifnoCol( 0 ) ;
+    ifno = ifnoCol( 0 ) ;
     os << LogIO::WARN
-       << "IFNO is not given. Using default IFNO: " << ifno_ << LogIO::POST ;
-  }
-  for ( uInt i = 0 ; i < nfile_ ; i++ ) {
-    Table taborg( infileList_[i] ) ;
-    TableExprNode node ;
-    if ( ifno != -1 || isMultiIF( taborg ) ) {
-      os << "apply selection on IFNO" << LogIO::POST ;
-      node = taborg.col("IFNO") == ifno_ ;
-    }
-    if ( scanlist_.size() > 0 ) {
-      os << "apply selection on SCANNO" << LogIO::POST ;
-      node = node && taborg.col("SCANNO").in( scanlist_ ) ;
-    }
-    if ( node.isNull() ) {
-      tableList_[i] = taborg ;
+       << "IFNO is not given. Using default IFNO: " << ifno << LogIO::POST ;
+  }
+//   tab = taborg( taborg.col("IFNO") == ifno ) ;
+  TableExprNode node ;
+  node = taborg.col("IFNO") == ifno ;
+  if ( scanlist_.size() > 0 ) {
+    node = node && taborg.col("SCANNO").in( scanlist_ ) ;
+  }
+  tab = taborg( node ) ;
+  if ( tab.nrow() == 0 ) {
+    LogIO os( LogOrigin("STGrid","selectData",WHERE) ) ;
+    os << LogIO::SEVERE
+       << "No corresponding rows for given selection: IFNO " << ifno 
+       << " SCANNO " << scanlist_ 
+       << LogIO::EXCEPTION ;
+  }
+}
+
+void STGrid::getData( Cube<Float> &spectra,
+                      Matrix<Double> &direction,
+                      Cube<uChar> &flagtra,
+                      Matrix<uInt> &rflag,
+                      Matrix<Float> &weight ) 
+{
+  Table tab ;
+  selectData( tab ) ;
+  updatePolList( tab ) ;
+//   cout << "npol_ = " << npol_ << endl ;
+//   cout << "nchan_ = " << nchan_ << endl ;
+//   cout << "nrow_ = " << nrow_ << endl ;
+  spectra.resize( npol_, nchan_, nrow_ ) ;
+  flagtra.resize( npol_, nchan_, nrow_ ) ;
+  rflag.resize( npol_, nrow_ ) ;
+  Cube<Float> tsys( npol_, nchan_, nrow_ ) ;
+  Matrix<Double> tint( npol_, nrow_ ) ;
+  // boolean for pointer access
+  Bool bsp, bfl, bfr, bts, bti ;
+  // pointer to the data
+  Float *sp_p = spectra.getStorage( bsp ) ;
+  uChar *fl_p = flagtra.getStorage( bfl ) ;
+  uInt *fr_p = rflag.getStorage( bfr ) ;
+  Float *ts_p = tsys.getStorage( bts ) ;
+  Double *ti_p = tint.getStorage( bti ) ;
+  // working pointer
+  Float *wsp_p = sp_p ;
+  uChar *wfl_p = fl_p ;
+  uInt *wfr_p = fr_p ;
+  Float *wts_p = ts_p ;
+  Double *wti_p = ti_p ;
+  uInt len = nchan_ * nrow_ ;
+  IPosition mshape( 2, nchan_, nrow_ ) ;
+  IPosition vshape( 1, nrow_ ) ;
+  for ( Int ipol = 0 ; ipol < npol_ ; ipol++ ) {
+    Table subt = tab( tab.col("POLNO") == pollist_[ipol] ) ;
+    ROArrayColumn<Float> spectraCol( subt, "SPECTRA" ) ;
+    ROArrayColumn<Double> directionCol( subt, "DIRECTION" ) ;
+    ROArrayColumn<uChar> flagtraCol( subt, "FLAGTRA" ) ;
+    ROScalarColumn<uInt> rflagCol( subt, "FLAGROW" ) ;
+    ROArrayColumn<Float> tsysCol( subt, "TSYS" ) ;
+    ROScalarColumn<Double> tintCol( subt, "INTERVAL" ) ;
+    Matrix<Float> spSlice( mshape, wsp_p, SHARE ) ; 
+    Matrix<uChar> flSlice( mshape, wfl_p, SHARE ) ;
+    Vector<uInt> frSlice( vshape, wfr_p, SHARE ) ;
+    spectraCol.getColumn( spSlice ) ;
+    flagtraCol.getColumn( flSlice ) ;
+    rflagCol.getColumn( frSlice ) ;
+    if ( ipol == 0 )
+      directionCol.getColumn( direction ) ;
+    Vector<Float> tmpF = tsysCol( 0 ) ;
+    Vector<Double> tmpD( vshape, wti_p, SHARE ) ;
+    Matrix<Float> tsSlice( mshape, wts_p, SHARE ) ;
+    if ( tmpF.nelements() == (uInt)nchan_ ) {
+      tsysCol.getColumn( tsSlice ) ;
     }
     else {
-      tableList_[i] = taborg( node ) ;
-    }
-    os << LogIO::DEBUGGING << "tableList_[" << i << "].nrow()=" << tableList_[i].nrow() << LogIO::POST ;
-    if ( tableList_[i].nrow() == 0 ) {
-      os << LogIO::SEVERE
-         << "No corresponding rows for given selection: IFNO " << ifno_ ;
-      if ( scanlist_.size() > 0 )
-        os << " SCANNO " << scanlist_ ;
-      os << LogIO::EXCEPTION ;
-    }
-  }
-}
-
-Bool STGrid::isMultiIF( Table &tab ) 
-{
-  ROScalarColumn<uInt> ifnoCol( tab, "IFNO" ) ;
-  Vector<uInt> ifnos = ifnoCol.getColumn() ;
-  return anyNE( ifnos, ifnos[0] ) ;
-}
-
-void STGrid::attach( Table &tab )
-{
-  // attach to table
-  spectraCol_.attach( tab, "SPECTRA" ) ;
-  flagtraCol_.attach( tab, "FLAGTRA" ) ;
-  directionCol_.attach( tab, "DIRECTION" ) ;
-  flagRowCol_.attach( tab, "FLAGROW" ) ;
-  tsysCol_.attach( tab, "TSYS" ) ;
-  intervalCol_.attach( tab, "INTERVAL" ) ;
-}
-
-Int STGrid::getDataChunk(
-                         IPosition const &wshape,
-                         IPosition const &vshape,
-                         IPosition const &dshape,
-                         Array<Complex> &spectra,
-                         Array<Double> &direction,
-                         Array<Int> &flagtra,
-                         Array<Int> &rflag,
-                         Array<Float> &weight ) 
-{
-  LogIO os( LogOrigin("STGrid","getDataChunk",WHERE) ) ;
-
-  Array<Float> spectraF_(wshape);
-  Array<uChar> flagtraUC_(wshape);
-  Array<uInt> rflagUI_(vshape);
-  Int nrow = getDataChunk( spectraF_, direction, flagtraUC_, rflagUI_, weight ) ;
-  if ( nrow < nchunk_ ) {
-    spectra.resize( spectraF_.shape() ) ;
-    flagtra.resize( flagtraUC_.shape() ) ;
-    rflag.resize( rflagUI_.shape() ) ;
-  }
-  double t0, t1 ;
-  t0 = mathutil::gettimeofday_sec() ;
-  convertArray( spectra, spectraF_ ) ;
-  toInt( flagtraUC_, flagtra ) ;
-  toInt( rflagUI_, rflag ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  eToInt = t1 - t0 ;
-  
-  return nrow ;
-}
-
-#if 0
-Int STGrid::getDataChunk( Array<Complex> &spectra,
-                          Array<Double> &direction,
-                          Array<Int> &flagtra,
-                          Array<Int> &rflag,
-                          Array<Float> &weight ) 
-{
-  LogIO os( LogOrigin("STGrid","getDataChunk",WHERE) ) ;
-  Int nrow = getDataChunk( spectraF_, direction, flagtraUC_, rflagUI_, weight ) ;
-  if ( nrow < nchunk_ ) {
-    spectra.resize( spectraF_.shape() ) ;
-    flagtra.resize( flagtraUC_.shape() ) ;
-    rflag.resize( rflagUI_.shape() ) ;
-  }
-  double t0, t1 ;
-  t0 = mathutil::gettimeofday_sec() ;
-  convertArray( spectra, spectraF_ ) ;
-  toInt( flagtraUC_, flagtra ) ;
-  toInt( rflagUI_, rflag ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  eToInt = t1 - t0 ;
-  
-  return nrow ;
-}
-#endif
-
-Int STGrid::getDataChunk( Array<Float> &spectra,
-                          Array<Double> &direction,
-                          Array<uChar> &flagtra,
-                          Array<uInt> &rflag,
-                          Array<Float> &weight ) 
-{
-  LogIO os( LogOrigin("STGrid","getDataChunk",WHERE) ) ;
-  Int nrow = spectra.shape()[1] ;
-  Int remainingRow = nrow_ - nprocessed_ ;
-  if ( remainingRow < nrow ) {
-    nrow = remainingRow ;
-    IPosition mshape( 2, nchan_, nrow ) ;
-    IPosition vshape( 1, nrow ) ;
-    spectra.resize( mshape ) ;
-    flagtra.resize( mshape ) ;
-    direction.resize( IPosition(2,2,nrow) ) ;
-    rflag.resize( vshape ) ;
-    weight.resize( mshape ) ;
-  }
-  // 2011/12/22 TN
-  // tsys shares its storage with weight
-  Array<Float> tsys( weight ) ;
-  Array<Double> tint( rflag.shape() ) ;
-
-  Vector<uInt> rflagVec( rflag ) ;
-  Vector<Double> tintVec( tint ) ;
-
-  RefRows rows( nprocessed_, nprocessed_+nrow-1, 1 ) ;
-  //os<<LogIO::DEBUGGING<<"nprocessed_="<<nprocessed_<<": rows.nrows()="<<rows.nrows()<<LogIO::POST ;
-  spectraCol_.getColumnCells( rows, spectra ) ;
-  flagtraCol_.getColumnCells( rows, flagtra ) ;
-  directionCol_.getColumnCells( rows, direction ) ;
-  flagRowCol_.getColumnCells( rows, rflagVec ) ;
-  intervalCol_.getColumnCells( rows, tintVec ) ;
-  Vector<Float> tsysTemp = tsysCol_( nprocessed_ ) ;
-  if ( tsysTemp.nelements() == (uInt)nchan_ )
-    tsysCol_.getColumnCells( rows, tsys ) ;
-  else
-    tsys = tsysTemp[0] ;
-
-  double t0,t1 ;
-  t0 = mathutil::gettimeofday_sec() ;
+      tsSlice = tmpF( 0 ) ;
+    }
+    tintCol.getColumn( tmpD ) ;
+
+    wsp_p += len ;
+    wfl_p += len ;
+    wfr_p += nrow_ ;
+    wts_p += len ;
+    wti_p += nrow_ ;
+  }
+  spectra.putStorage( sp_p, bsp ) ;
+  flagtra.putStorage( fl_p, bfl ) ;
+  rflag.putStorage( fr_p, bfr ) ;
+  tsys.putStorage( ts_p, bts ) ;
+  tint.putStorage( ti_p, bti ) ;
+
   getWeight( weight, tsys, tint ) ;
-  t1 = mathutil::gettimeofday_sec() ;
-  eGetWeight += t1-t0 ;
-
-  nprocessed_ += nrow ;
-  
-  return nrow ;
-}
-
-void STGrid::setupArray() 
-{
-  LogIO os( LogOrigin("STGrid","setupArray",WHERE) ) ;
-  ROScalarColumn<uInt> polnoCol( tableList_[0], "POLNO" ) ;
+}
+
+void STGrid::updatePolList( Table &tab ) 
+{
+  ROScalarColumn<uInt> polnoCol( tab, "POLNO" ) ;
   Vector<uInt> pols = polnoCol.getColumn() ;
-  //os << pols << LogIO::POST ;
   Vector<uInt> pollistOrg ;
-  npolOrg_ = 0 ;
+  uInt npolOrg = 0 ;
   uInt polno ;
   for ( uInt i = 0 ; i < polnoCol.nrow() ; i++ ) {
@@ -1365 +544 @@
     polno = pols( i ) ; 
     if ( allNE( pollistOrg, polno ) ) {
-      pollistOrg.resize( npolOrg_+1, True ) ;
-      pollistOrg[npolOrg_] = polno ;
-      npolOrg_++ ;
+      pollistOrg.resize( npolOrg+1, True ) ;
+      pollistOrg[npolOrg] = polno ;
+      npolOrg++ ;
     }
   }
@@ -1386 +565 @@
   npol_ = pollist_.size() ;
   if ( npol_ == 0 ) {
+    LogIO os( LogOrigin("STGrid","updatePolList",WHERE) ) ;
     os << LogIO::SEVERE << "Empty pollist" << LogIO::EXCEPTION ;
   }
-  rows_.resize( nfile_ ) ;
-  for ( uInt i = 0 ; i < nfile_ ; i++ ) {
-    rows_[i] = tableList_[i].nrow() / npolOrg_ ;
-    //if ( nrow_ < rows_[i] ) 
-    //  nrow_ = rows_[i] ;
-  }
-  flagtraCol_.attach( tableList_[0], "FLAGTRA" ) ;
-  nchan_ = flagtraCol_( 0 ).nelements() ;
+  nrow_ = tab.nrow() / npolOrg ;
+  ROArrayColumn<uChar> tmpCol( tab, "FLAGTRA" ) ;
+  nchan_ = tmpCol( 0 ).nelements() ;
+//   LogIO os( LogOrigin("STGrid","updatePolList",WHERE) ) ;
 //   os << "npol_ = " << npol_ << "(" << pollist_ << ")" << endl 
 //      << "nchan_ = " << nchan_ << endl 
@@ -1401 +577 @@
 }
 
-void STGrid::getWeight( Array<Float> &w,
-                              Array<Float> &tsys,
-                              Array<Double> &tint ) 
+void STGrid::getWeight( Matrix<Float> &w,
+                        Cube<Float> &tsys,
+                        Matrix<Double> &tint ) 
 {
   LogIO os( LogOrigin("STGrid","getWeight",WHERE) ) ;
-
-  // 2011/12/22 TN
-  // w (weight) and tsys share storage
-  IPosition refShape = tsys.shape() ;
-  Int nchan = refShape[0] ;
-  Int nrow = refShape[1] ;
-//   os << "nchan=" << nchan << ", nrow=" << nrow << LogIO::POST ;
-//   os << "w.shape()=" << w.shape() << endl
-//      << "tsys.shape()=" << tsys.shape() << endl
-//      << "tint.shape()=" << tint.shape() << LogIO::POST ;
+  double t0, t1 ;
+  t0 = mathutil::gettimeofday_sec() ;
+  // resize
+  w.resize( nchan_, nrow_ ) ;
 
   // set weight
+  Bool warn = False ;
   if ( wtype_.compare( "UNIFORM" ) == 0 ) {
     w = 1.0 ;
   }
   else if ( wtype_.compare( "TINT" ) == 0 ) {
+    if ( npol_ > 1 ) warn = True ;
     Bool b0, b1 ;
     Float *w_p = w.getStorage( b0 ) ;
@@ -1427 +599 @@
     const Double *ti_p = tint.getStorage( b1 ) ;
     const Double *w1_p = ti_p ;
-    for ( Int irow = 0 ; irow < nrow ; irow++ ) {
-      for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
-        *w0_p = *w1_p ;
+    for ( Int irow = 0 ; irow < nrow_ ; irow++ ) {
+      Float val = (Float)(polMean( w1_p )) ;
+      for ( Int ichan = 0 ; ichan < nchan_ ; ichan++ ) {
+        *w0_p = val ;
         w0_p++ ;
       }
-      w1_p++ ;
     }
     w.putStorage( w_p, b0 ) ;
@@ -1438 +610 @@
   }
   else if ( wtype_.compare( "TSYS" ) == 0 ) {
-    Bool b0 ;
+    if ( npol_ > 1 ) warn = True ;
+    Bool b0, b1 ;
     Float *w_p = w.getStorage( b0 ) ;
     Float *w0_p = w_p ;
-    for ( Int irow = 0 ; irow < nrow ; irow++ ) {
-      for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
-        Float temp = *w0_p ;
-        *w0_p = 1.0 / ( temp * temp ) ;
+    const Float *ts_p = tsys.getStorage( b1 ) ;
+    const Float *w1_p = ts_p ;
+    for ( Int irow = 0 ; irow < nrow_ ; irow++ ) {
+      for ( Int ichan = 0 ; ichan < nchan_ ; ichan++ ) {
+        Float val = polMean( w1_p ) ;
+        *w0_p = 1.0 / ( val * val ) ;
         w0_p++ ;
       }
     }
     w.putStorage( w_p, b0 ) ;
+    tsys.freeStorage( ts_p, b1 ) ;
   }
   else if ( wtype_.compare( "TINTSYS" ) == 0 ) {
-    Bool b0, b1 ;
+    if ( npol_ > 1 ) warn = True ;
+    Bool b0, b1, b2 ;
     Float *w_p = w.getStorage( b0 ) ;
     Float *w0_p = w_p ;
     const Double *ti_p = tint.getStorage( b1 ) ;
     const Double *w1_p = ti_p ;
-    for ( Int irow = 0 ; irow < nrow ; irow++ ) {
-      Float interval = *w1_p ;
-      for ( Int ichan = 0 ; ichan < nchan ; ichan++ ) {
-        Float temp = *w0_p ;
+    const Float *ts_p = tsys.getStorage( b2 ) ;
+    const Float *w2_p = ts_p ;
+    for ( Int irow = 0 ; irow < nrow_ ; irow++ ) {
+      Float interval = (Float)(polMean( w1_p )) ;
+      for ( Int ichan = 0 ; ichan < nchan_ ; ichan++ ) {
+        Float temp = polMean( w2_p ) ;
         *w0_p = interval / ( temp * temp ) ;
         w0_p++ ;
       }
-      w1_p++ ;
     }
     w.putStorage( w_p, b0 ) ;
     tint.freeStorage( ti_p, b1 ) ;
+    tsys.freeStorage( ts_p, b2 ) ;
   }
   else {
     //LogIO os( LogOrigin("STGrid", "getWeight", WHERE) ) ;
-    //os << LogIO::WARN << "Unsupported weight type '" << wtype_ << "', apply UNIFORM weight" << LogIO::POST ;
+    os << LogIO::WARN << "Unsupported weight type '" << wtype_ << "', apply UNIFORM weight" << LogIO::POST ;
     w = 1.0 ;
   }
-}
-
-void STGrid::toInt( Array<uChar> &u, Array<Int> &v ) 
-{
-  uInt len = u.nelements() ;
+
+  if ( npol_ > 1 ) {
+    //LogIO os( LogOrigin("STGrid", "getWeight", WHERE) ) ;
+    os << LogIO::WARN << "STGrid doesn't support assigning polarization-dependent weight. Use averaged weight over polarization." << LogIO::POST ;
+  }
+  t1 = mathutil::gettimeofday_sec() ;
+  os << "getWeight: elapsed time is " << t1-t0 << " sec" << LogIO::POST ;
+}
+
+Float STGrid::polMean( const Float *p ) 
+{
+  Float v = 0.0 ;
+  for ( Int i = 0 ; i < npol_ ; i++ ) {
+    v += *p ;
+    p++ ;
+  }
+  v /= npol_ ;
+  return v ;
+}
+
+Double STGrid::polMean( const Double *p ) 
+{
+  Double v = 0.0 ;
+  for ( Int i = 0 ; i < npol_ ; i++ ) {
+    v += *p ;
+    p++ ;
+  }
+  v /= npol_ ;
+  return v ;
+}
+
+void STGrid::toInt( Array<uChar> *u, Array<Int> *v ) 
+{
+  uInt len = u->nelements() ;
   Int *int_p = new Int[len] ;
   Bool deleteIt ;
-  const uChar *data_p = u.getStorage( deleteIt ) ;
+  const uChar *data_p = u->getStorage( deleteIt ) ;
   Int *i_p = int_p ;
   const uChar *u_p = data_p ;
@@ -1488 +696 @@
     u_p++ ;
   }
-  u.freeStorage( data_p, deleteIt ) ;
-  v.takeStorage( u.shape(), int_p, TAKE_OVER ) ;
-}
-
-void STGrid::toInt( Array<uInt> &u, Array<Int> &v ) 
-{
-  uInt len = u.nelements() ;
+  u->freeStorage( data_p, deleteIt ) ;
+  v->takeStorage( u->shape(), int_p, TAKE_OVER ) ;
+}
+
+void STGrid::toInt( Array<uInt> *u, Array<Int> *v ) 
+{
+  uInt len = u->nelements() ;
   Int *int_p = new Int[len] ;
   Bool deleteIt ;
-  const uInt *data_p = u.getStorage( deleteIt ) ;
+  const uInt *data_p = u->getStorage( deleteIt ) ;
   Int *i_p = int_p ;
   const uInt *u_p = data_p ;
@@ -1505 +713 @@
     u_p++ ;
   }
-  u.freeStorage( data_p, deleteIt ) ;
-  v.takeStorage( u.shape(), int_p, TAKE_OVER ) ;
-}
-
-void STGrid::toPixel( Array<Double> &world, Array<Double> &pixel )
+  u->freeStorage( data_p, deleteIt ) ;
+  v->takeStorage( u->shape(), int_p, TAKE_OVER ) ;
+}
+
+void STGrid::toPixel( Matrix<Double> &world, Matrix<Double> &pixel )
 {
   // gridding will be done on (nx_+2*convSupport_) x (ny_+2*convSupport_) 
   // grid plane to avoid unexpected behavior on grid edge
-  Block<Double> pixc( 2 ) ;
-  pixc[0] = Double( nx_-1 ) * 0.5 ;
-  pixc[1] = Double( ny_-1 ) * 0.5 ;
-//   pixc[0] = Double( nx_+2*convSupport_-1 ) * 0.5 ;
-//   pixc[1] = Double( ny_+2*convSupport_-1 ) * 0.5 ;
+  Vector<Double> pixc( 2 ) ;
+  pixc(0) = Double( nx_-1 ) * 0.5 ;
+  pixc(1) = Double( ny_-1 ) * 0.5 ;
+//   pixc(0) = Double( nx_+2*convSupport_-1 ) * 0.5 ;
+//   pixc(1) = Double( ny_+2*convSupport_-1 ) * 0.5 ;
   uInt nrow = world.shape()[1] ;
-  Bool bw, bp ;
-  const Double *w_p = world.getStorage( bw ) ;
-  Double *p_p = pixel.getStorage( bp ) ;
-  const Double *ww_p = w_p ;
-  Double *wp_p = p_p ;
-  for ( uInt i = 0 ; i < nrow ; i++ ) {
-    *wp_p = pixc[0] + ( *ww_p - center_[0] ) / cellx_ ;
-    wp_p++ ;
-    ww_p++ ;
-    *wp_p = pixc[1] + ( *ww_p - center_[1] ) / celly_ ;
-    wp_p++ ;
-    ww_p++ ;
-  }
-  world.freeStorage( w_p, bw ) ;
-  pixel.putStorage( p_p, bp ) ;  
+  Vector<Double> cell( 2 ) ;
+  cell(0) = cellx_ ;
+  cell(1) = celly_ ;
+  for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
+    for ( uInt i = 0 ; i < 2 ; i++ ) {
+      pixel( i, irow ) = pixc(i) + ( world(i, irow) - center_(i) ) / cell(i) ;
+    }
+  }
+//   String gridfile = "grid."+convType_+"."+String::toString(convSupport_)+".dat" ; 
+//   ofstream ofs( gridfile.c_str(), ios::out ) ; 
+//   ofs << "center " << center_(0) << " " << pixc(0) 
+//       << " " << center_(1) << " " << pixc(1) << endl ;
+//   for ( uInt irow = 0 ; irow < nrow ; irow++ ) {
+//     ofs << irow ;
+//     for ( uInt i = 0 ; i < 2 ; i++ ) {
+//       ofs << " " << world(i, irow) << " " << pixel(i, irow) ;
+//     }
+//     ofs << endl ;
+//   }
+//   ofs.close() ;
 }
 
@@ -1599 +812 @@
     // to take into account Gaussian tail
     if ( convSupport_ < 0 )
-      convSupport_ = 4 ; // 1 * 4
+      convSupport_ = 12 ; // 3 * 4
     else {
       convSupport_ = userSupport_ * 4 ;
@@ -1626 +839 @@
   String outfile_ ;
   if ( outfile.size() == 0 ) {
-    if ( infileList_[0].lastchar() == '/' ) {
-      outfile_ = infileList_[0].substr( 0, infileList_[0].size()-1 ) ;
+    if ( infile_.lastchar() == '/' ) {
+      outfile_ = infile_.substr( 0, infile_.size()-1 ) ;
     }
     else {
-      outfile_ = infileList_[0] ;
+      outfile_ = infile_ ;
     }
     outfile_ += ".grid" ;
@@ -1684 +897 @@
 
   t1 = mathutil::gettimeofday_sec() ;
-  os << LogIO::DEBUGGING << "saveData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
-
-  fillMainColumns( tab ) ;
-
+  os << "saveData: elapsed time is " << t1-t0 << " sec." << LogIO::POST ; 
+  
   return outfile_ ;
 }
@@ -1693 +904 @@
 void STGrid::prepareTable( Table &tab, String &name ) 
 {
-  Table t( infileList_[0], Table::Old ) ;
+  Table t( infile_, Table::Old ) ;
   t.deepCopy( name, Table::New, False, t.endianFormat(), True ) ;
   tab = Table( name, Table::Update ) ;
-  // 2012/02/13 TN
-  // explicitly copy subtables since no rows including subtables are 
-  // copied by Table::deepCopy with noRows=True
-  TableCopy::copySubTables( tab, t ) ;
-}
-
-void STGrid::fillMainColumns( Table &tab ) 
-{
-  // values for fill
-  Table t( infileList_[0], Table::Old ) ;
-  Table tsel = t( t.col( "IFNO" ) == (uInt)ifno_, 1 ) ;
-  ROTableRow row( tsel ) ;
-  row.get( 0 ) ;
-  const TableRecord &rec = row.record() ;
-  uInt freqId = rec.asuInt( "FREQ_ID" ) ;
-  uInt molId = rec.asuInt( "MOLECULE_ID" ) ;
-  uInt tcalId = rec.asuInt( "TCAL_ID" ) ;
-  uInt focusId = rec.asuInt( "FOCUS_ID" ) ;
-  uInt weatherId = rec.asuInt( "WEATHER_ID" ) ;
-  String srcname = rec.asString( "SRCNAME" ) ;
-  String fieldname = rec.asString( "FIELDNAME" ) ;
-  Vector<Float> defaultTsys( 1, 1.0 ) ;
-  // @todo how to set flagtra for gridded spectra?
-  Vector<uChar> flagtra = rec.asArrayuChar( "FLAGTRA" ) ;
-  flagtra = (uChar)0 ;
-  Float opacity = rec.asFloat( "OPACITY" ) ;
-  Double srcvel = rec.asDouble( "SRCVELOCITY" ) ;
-  Vector<Double> srcpm = rec.asArrayDouble( "SRCPROPERMOTION" ) ;
-  Vector<Double> srcdir = rec.asArrayDouble( "SRCDIRECTION" ) ;
-  Vector<Double> scanrate = rec.asArrayDouble( "SCANRATE" ) ;
-  Double time = rec.asDouble( "TIME" ) ;
-  Double interval = rec.asDouble( "INTERVAL" ) ;
-
-  // fill columns
-  Int nrow = tab.nrow() ;
-  ScalarColumn<uInt> scannoCol( tab, "SCANNO" ) ;
-  ScalarColumn<uInt> ifnoCol( tab, "IFNO" ) ;
-  ScalarColumn<uInt> freqIdCol( tab, "FREQ_ID" ) ;
-  ScalarColumn<uInt> molIdCol( tab, "MOLECULE_ID" ) ;
-  ScalarColumn<uInt> tcalidCol( tab, "TCAL_ID" ) ;
-  ScalarColumn<Int> fitidCol( tab, "FIT_ID" ) ;
-  ScalarColumn<uInt> focusidCol( tab, "FOCUS_ID" ) ;
-  ScalarColumn<uInt> weatheridCol( tab, "WEATHER_ID" ) ;
-  ArrayColumn<uChar> flagtraCol( tab, "FLAGTRA" ) ;
-  ScalarColumn<uInt> rflagCol( tab, "FLAGROW" ) ;
-  ArrayColumn<Float> tsysCol( tab, "TSYS" ) ;
-  ScalarColumn<String> srcnameCol( tab, "SRCNAME" ) ;
-  ScalarColumn<String> fieldnameCol( tab, "FIELDNAME" ) ;
-  ScalarColumn<Int> srctypeCol( tab, "SRCTYPE" ) ;
-  ScalarColumn<Float> opacityCol( tab, "OPACITY" ) ;
-  ScalarColumn<Double> srcvelCol( tab, "SRCVELOCITY" ) ;
-  ArrayColumn<Double> srcpmCol( tab, "SRCPROPERMOTION" ) ;
-  ArrayColumn<Double> srcdirCol( tab, "SRCDIRECTION" ) ;
-  ArrayColumn<Double> scanrateCol( tab, "SCANRATE" ) ;
-  ScalarColumn<Double> timeCol( tab, "TIME" ) ;
-  ScalarColumn<Double> intervalCol( tab, "INTERVAL" ) ;
-  for ( Int i = 0 ; i < nrow ; i++ ) {
-    scannoCol.put( i, (uInt)i ) ;
-    ifnoCol.put( i, (uInt)ifno_ ) ;
-    freqIdCol.put( i, freqId ) ;
-    molIdCol.put( i, molId ) ;
-    tcalidCol.put( i, tcalId ) ;
-    fitidCol.put( i, -1 ) ;
-    focusidCol.put( i, focusId ) ;
-    weatheridCol.put( i, weatherId ) ;
-    flagtraCol.put( i, flagtra ) ;
-    rflagCol.put( i, 0 ) ;
-    tsysCol.put( i, defaultTsys ) ;
-    srcnameCol.put( i, srcname ) ;
-    fieldnameCol.put( i, fieldname ) ;
-    srctypeCol.put( i, (Int)SrcType::PSON ) ;
-    opacityCol.put( i, opacity ) ;
-    srcvelCol.put( i, srcvel ) ;
-    srcpmCol.put( i, srcpm ) ;
-    srcdirCol.put( i, srcdir ) ;
-    scanrateCol.put( i, scanrate ) ;
-    timeCol.put( i, time ) ;
-    intervalCol.put( i, interval ) ;
-  }
-}
-
-}
+}
+}