Changeset 2126

Timestamp:

04/08/11 21:06:20 (13 years ago)

Author:

WataruKawasaki

Message:

merged bug fixes from trunk (r2124)

Location:

branches/casa-prerelease/pre-asap

Files:

• . (modified) (1 prop)
• Makefile (modified) (1 prop)
• SConstruct (modified) (1 prop)
• apps (modified) (1 prop)
• external-alma (modified) (1 prop)
• external-alma/atnf/pks/pks_maths.cc (modified) (1 prop)
• getsvnrev.sh (modified) (1 prop)
• python (modified) (1 prop)
• src (modified) (1 prop)
• src/MathUtils.cpp (modified) (3 diffs)
• src/RowAccumulator.cpp (modified) (8 diffs)
• src/RowAccumulator.h (modified) (4 diffs)
• src/SConscript (modified) (1 prop)
• src/STMath.cpp (modified) (8 diffs)
• src/python_Scantable.cpp (modified) (1 diff)

branches/casa-prerelease/pre-asap/src/MathUtils.cpp

@@ -30 +30 @@
 //#---------------------------------------------------------------------------
 
+#include <iostream>
 #include <casa/aips.h>
 #include <casa/Arrays/Vector.h>
@@ -161 +162 @@
   out.resize(in.nelements());
   out[0] = in[0];out[out.nelements()-1] = in[in.nelements()-1];
+  //--s
+  cout << "out[0]=" << out[0] << ";" << endl;
+  //--e
   outmask.resize(mask.nelements());
-  outmask = False;
+  outmask[0] = mask[0]; outmask[outmask.nelements()-1] = mask[mask.nelements()-1];
   uInt m;Vector<Float>* w;
   for (uInt i=1; i < out.nelements()-1;++i) {
@@ -169 +173 @@
     if (weighted[m]) {
       out[i] = (*w)[0]*in[i-1] + (*w)[1]*in[i] + (*w)[2]*in[i+1];
-      outmask[i] = True;
+      //outmask[i] = True;
     } else { // mask it
       out[i] = in[i];//use arbitrary value
-      outmask[i] = False;
+      //outmask[i] = False;
     }
+    //--s
+    if (i < 13) {
+      cout << "mask[" << i-1 << "-" << i+1 << "]=[" << mask[i-1] << ", " << mask[i] << ", " << mask[i+1] << "];  ";
+      cout << "m="   << m   << ";  ";
+      cout << "weights[" << m << "]="   << weights[m]   << ";  ";
+      cout << "in[" << i-1 << "-" << i+1 << "]=[" << in[i-1] << ", " << in[i] << ", " << in[i+1] << "] --> ";
+      cout << out[i]  << ";" << endl;
+      cout << "-------------------------------" << endl;
+    }
+    //--e
+    outmask[i] = mask[i];
   }
 }

branches/casa-prerelease/pre-asap/src/RowAccumulator.cpp

@@ -10 +10 @@
 //
 //
+
+#include <iostream>
+
 #include <casa/iomanip.h>
-
 #include <casa/Arrays/MaskArrMath.h>
 #include <casa/Arrays/MaskArrLogi.h>
@@ -18 +20 @@
 #include "RowAccumulator.h"
 
-
 using namespace casa;
 using namespace asap;
 
-RowAccumulator::RowAccumulator(WeightType wt) :
-  weightType_(wt),
-  initialized_(False)
+RowAccumulator::RowAccumulator(WeightType wt) : weightType_(wt), initialized_(False)
 {
   reset();
@@ -34 +33 @@
 
 
-void RowAccumulator::add( const Vector< Float >& v,
-                          const Vector< Bool >& m,
-                          const Vector< Float >& tsys,
-                          Double interval,
-                          Double time )
-{
-  if (!initialized_) {
-    Vector<Float> dummy(v.nelements(), 0.0);
-    Vector<Bool> dummymsk(m.nelements(), True);
-    spectrum_.setData(dummy, dummymsk);
-    n_.setData(Vector<uInt>(v.nelements(), 0), dummymsk);
-    weightSum_.setData(Vector<Float>(v.nelements(), 0.0), dummymsk);
-    tsysSum_.resize(tsys.nelements()); tsysSum_=0.0;
-  }
-  // add spectrum related weights, so far it is variance only.
-  Float totalweight = 1.0;
-
-  // only add these if not everything masked
-  if ( !allEQ(m, False) ) {
-    totalweight *= addTsys(tsys);
-    totalweight *= addInterval(interval);
-    addTime(time);
-  }
-  addSpectrum(v, m, totalweight);
-  initialized_ = True;
-}
-
-void RowAccumulator::addSpectrum( const Vector< Float >& v,
-                                  const Vector< Bool >& m,
-                                  Float weight)
-{
-  Float totalweight = weight;
-  MaskedArray<Float> data(v,m);
-  if ( weightType_ == asap::W_VAR ) {
-    if (m.nelements() == userMask_.nelements()) {
-      Float fac = 1.0/variance(data(userMask_));
-      totalweight *= fac;
-    } else {
-      Float fac = 1.0/variance(data);
-      totalweight *= fac;
-    }
-  }
-  data *= totalweight;
-  MaskedArray<Float> wadd(Vector<Float>(m.nelements(),totalweight), m);
-  weightSum_ += wadd;
-  spectrum_ += data;
-  const MaskedArray<uInt> inc(Vector<uInt>(m.nelements(),1), m);
-  n_ += inc;
-}
-
-Float RowAccumulator::addTsys( const casa::Vector< casa::Float > & v )
+void RowAccumulator::reset(const uInt size, const uInt tsysSize)
+{
+  Vector<Bool> m(size, True);
+
+  spectrum_.setData(Vector<Float>(size, 0.0), Vector<Bool>(size, True));
+  spectrumNoMask_.setData(Vector<Float>(size, 0.0), Vector<Bool>(size, True));
+
+  n_.setData(Vector<uInt>(size, 0), Vector<Bool>(size, True));
+  nNoMask_.setData(Vector<uInt>(size, 0), Vector<Bool>(size, True));
+
+  weightSum_.setData(Vector<Float>(size, 0.0), Vector<Bool>(size, True));
+  weightSumNoMask_.setData(Vector<Float>(size, 0.0), Vector<Bool>(size, True));
+
+  tsysSum_.resize(tsysSize); tsysSum_=0.0;
+  tsysSumNoMask_.resize(tsysSize); tsysSumNoMask_=0.0;
+
+  intervalSum_ = 0.0;
+  intervalSumNoMask_ = 0.0;
+
+  timeSum_ = 0.0;
+  timeSumNoMask_ = 0.0;
+}
+
+void RowAccumulator::initialize(const uInt size, const uInt tsysSize)
+{
+    reset(size, tsysSize);
+    initialized_ = True;
+}
+
+void RowAccumulator::add(const Vector<Float>& v,
+                         const Vector<Bool>& m,
+                         const Vector<Float>& tsys,
+                         const Double interval,
+                         const Double time)
+{
+  uInt size = v.nelements();
+  //if (size != m.nelements()) raiseError;
+  if (!initialized_) initialize(size, tsys.nelements());
+
+  addSpectrum(v, m, tsys, interval, time);
+}
+
+void RowAccumulator::addSpectrum(const Vector<Float>& v,
+                                 const Vector<Bool>& m,
+                                 const Vector<Float>& tsys,
+                                 const Double interval,
+                                 const Double time)
+{
+  doAddSpectrum(v, m, tsys, interval, time, False);
+  doAddSpectrum(v, m, tsys, interval, time, True);  // CAS-2776
+}
+
+void RowAccumulator::doAddSpectrum(const Vector<Float>& v,
+                                   const Vector<Bool>& m,
+                                   const Vector<Float>& tsys,
+                                   const Double interval,
+                                   const Double time,
+                                   const Bool ignoreMask)
+{
+  Vector<Float> vUse = v.copy();
+  Vector<Bool> mUse = m.copy();
+  if (ignoreMask) mUse = !mUse;
+
+  MaskedArray<Float> vadd(vUse, mUse);
+  Float totalWeight = getTotalWeight(vadd, tsys, interval, time, ignoreMask);
+  vadd *= totalWeight;
+  MaskedArray<Float> wadd(Vector<Float>(mUse.nelements(), totalWeight), mUse);
+  MaskedArray<uInt> inc(Vector<uInt>(mUse.nelements(), 1), mUse);
+
+  if (ignoreMask) {
+    spectrumNoMask_ += vadd;
+    weightSumNoMask_ += wadd;
+    nNoMask_ += inc;
+  } else {
+    spectrum_ += vadd;
+    weightSum_ += wadd;
+    n_ += inc;
+  }
+
+  //
+  cout << "***" << endl;
+  Vector<Float> spe = spectrum_.getArray();
+  Vector<Float> spe0 = spectrumNoMask_.getArray();
+  Vector<Float> wei = weightSum_.getArray();
+  Vector<Float> wei0 = weightSumNoMask_.getArray();
+  Vector<uInt>  n = n_.getArray();
+  Vector<uInt>  n0 = nNoMask_.getArray();
+  cout << "S__" << "[" << spe[0] << "][" << spe[1] << "][" << spe[2] << "][" << spe[3] << "][" << spe[4] << "][" << spe[5] << "][" << spe[6] << "][" << spe[7] << "][" << spe[8] << "][" << spe[9] << "][" << spe[10] << "][" << spe[11] << "][" << spe[12] << "]" << endl;
+  cout << "S0_" << "[" << spe0[0] << "][" << spe0[1] << "][" << spe0[2] << "][" << spe0[3] << "][" << spe0[4] << "][" << spe0[5] << "][" << spe0[6] << "][" << spe0[7] << "][" << spe0[8] << "][" << spe0[9] << "][" << spe0[10] << "][" << spe0[11] << "][" << spe0[12] << "]" << endl;
+  cout << "W__" << "[" << wei[0] << "][" << wei[1] << "][" << wei[2] << "][" << wei[3] << "][" << wei[4] << "][" << wei[5] << "][" << wei[6] << "][" << wei[7] << "][" << wei[8] << "][" << wei[9] << "][" << wei[10] << "][" << wei[11] << "][" << wei[12] << "]" << endl;
+  cout << "W0_" << "[" << wei0[0] << "][" << wei0[1] << "][" << wei0[2] << "][" << wei0[3] << "][" << wei0[4] << "][" << wei0[5] << "][" << wei0[6] << "][" << wei0[7] << "][" << wei0[8] << "][" << wei0[9] << "][" << wei0[10] << "][" << wei0[11] << "][" << wei0[12] << "]" << endl;
+  cout << "N__" << "[" << n[0] << "][" << n[1] << "][" << n[2] << "][" << n[3] << "][" << n[4] << "][" << n[5] << "][" << n[6] << "][" << n[7] << "][" << n[8] << "][" << n[9] << "][" << n[10] << "][" << n[11] << "][" << n[12] << "]" << endl;
+  cout << "N0_" << "[" << n0[0] << "][" << n0[1] << "][" << n0[2] << "][" << n0[3] << "][" << n0[4] << "][" << n0[5] << "][" << n0[6] << "][" << n0[7] << "][" << n0[8] << "][" << n0[9] << "][" << n0[10] << "][" << n0[11] << "][" << n0[12] << "]" << endl;
+  cout << "***" << endl;
+  //
+}
+
+Float RowAccumulator::getTotalWeight(const MaskedArray<Float>& data,
+                                     const Vector<Float>& tsys,
+                                     const Double interval,
+                                     const Double time,
+                                     const Bool ignoreMask)
+{
+  Float totalWeight = 1.0;
+
+  Vector<Bool> m = data.getMask();
+  if (!allEQ(m, False)) {  // only add these if not everything masked
+    totalWeight *= addTsys(tsys, ignoreMask);
+    totalWeight *= addInterval(interval, ignoreMask);
+    addTime(time, ignoreMask);
+  }
+
+  if (weightType_ == W_VAR) {
+    Float fac = 1.0/variance(data);
+    if (!ignoreMask && (m.nelements() == userMask_.nelements()))
+      fac = 1.0/variance(data(userMask_));
+
+    totalWeight *= fac;
+  }
+
+  return totalWeight;
+}
+
+Float RowAccumulator::addTsys(const Vector<Float>& v, Bool ignoreMask)
 {
   // @fixme this assume tsys is the same for all channels
 
   Float w = 1.0;
-  tsysSum_ += v[0];
-  if ( weightType_ == asap::W_TSYS  || weightType_ == asap::W_TINTSYS ) {
+  if (ignoreMask) {
+    tsysSumNoMask_ += v[0];
+  } else {
+    tsysSum_ += v[0];
+  }
+  if ( weightType_ == W_TSYS  || weightType_ == W_TINTSYS ) {
     w /= (v[0]*v[0]);
   }
@@ -96 +172 @@
 }
 
-void asap::RowAccumulator::addTime( casa::Double t )
-{
-  timeSum_ += t;
-}
-
-Float asap::RowAccumulator::addInterval( casa::Double inter )
+void RowAccumulator::addTime(Double t, Bool ignoreMask)
+{
+  if (ignoreMask) {
+    timeSumNoMask_ += t;
+  } else {
+    timeSum_ += t;
+  }
+}
+
+Float RowAccumulator::addInterval(Double inter, Bool ignoreMask)
 {
   Float w = 1.0;
-  intervalSum_ += inter;
-  if ( weightType_ == asap::W_TINT || weightType_ == asap::W_TINTSYS ) {
+  if (ignoreMask) {
+    intervalSumNoMask_ += inter;
+  } else {
+    intervalSum_ += inter;
+  }
+  if (weightType_ == W_TINT || weightType_ == W_TINTSYS) {
     w /= Float(inter);
   }
@@ -111 +195 @@
 }
 
-void asap::RowAccumulator::reset( )
-{
-  initialized_ = False;
-  intervalSum_ = 0.0;
-  tsysSum_.resize();
-  timeSum_ = 0.0;
-}
-
-casa::Vector< casa::Float > RowAccumulator::getSpectrum( ) const
+Vector<Float> RowAccumulator::getSpectrum() const
 {
   return (spectrum_/weightSum_).getArray();
 }
 
-casa::Double asap::RowAccumulator::getTime( ) const
-{
-  uInt n = max(n_);
-  return timeSum_/Float(n);
-}
-
-casa::Double asap::RowAccumulator::getInterval( ) const
+Double RowAccumulator::getTime() const
+{
+  return timeSum_/Float(max(n_));
+}
+
+Double RowAccumulator::getInterval() const
 {
   return intervalSum_;
 }
 
-casa::Vector< casa::Bool > RowAccumulator::getMask( ) const
+Vector<Bool> RowAccumulator::getMask() const
 {
   // Return the "total" mask - False where no points have been accumulated.
@@ -141 +216 @@
 }
 
-casa::Vector< casa::Float > asap::RowAccumulator::getTsys( ) const
+Vector<Float> RowAccumulator::getTsys() const
 {
   // @fixme this assumes tsys.nelements() == 1
@@ -147 +222 @@
 }
 
-void asap::RowAccumulator::setUserMask( const casa::Vector< casa::Bool > & m )
+void RowAccumulator::setUserMask(const Vector<Bool>& m)
 {
   userMask_.resize();
@@ -154 +229 @@
 
 // Added by TT  check the state of RowAccumulator
-casa::Bool RowAccumulator::state() const
+Bool RowAccumulator::state() const
 {
   return initialized_;
 }
 
+void RowAccumulator::replaceNaN()
+{
+  Vector<Float> v = spectrum_.getArray();
+  Vector<Float> w = weightSum_.getArray();
+  Vector<Float> vRef = spectrumNoMask_.getArray();
+  Vector<Float> wRef = weightSumNoMask_.getArray();
+
+  //-------
+  cout << "SB-";
+  for (uInt i=0; i<13; ++i) {
+    cout << "[" << v[i] << "]";
+  }
+  cout << endl;
+  cout << "WB-";
+  for (uInt i=0; i<13; ++i) {
+    cout << "[" << w[i] << "]";
+  }
+  cout << endl;
+  //-------
+
+  for (uInt i = 0; i < v.size(); ++i) {
+    if (w[i] == 0.0) {
+      v[i] = vRef[i];
+      w[i] = wRef[i];
+    }
+  }
+
+  spectrum_.setData(v, Vector<Bool>(v.nelements(), True));
+  weightSum_.setData(w, Vector<Bool>(w.nelements(), True));
+
+  //-------
+  cout << "SA-";
+  for (uInt i=0; i<13; ++i) {
+    cout << "[" << v[i] << "]";
+  }
+  cout << endl;
+  cout << "WA-";
+  for (uInt i=0; i<13; ++i) {
+    cout << "[" << w[i] << "]";
+  }
+  cout << endl;
+  //-------
+}

branches/casa-prerelease/pre-asap/src/RowAccumulator.h

@@ -13 +13 @@
 #define ASAPROWACCUMULATOR_H
 
+#include <math.h>
 #include <casa/aips.h>
 #include <casa/Arrays/Vector.h>
@@ -48 +49 @@
            const casa::Vector<casa::Bool>& m,
            const casa::Vector<casa::Float>& tsys,
-           casa::Double interval,
-           casa::Double time);
+           const casa::Double interval,
+           const casa::Double time);
   /**
     * Also set a user mask which get combined with the individual masks
@@ -84 +85 @@
     * Reset the acummulator to the state at construction.
     */
-  void reset();
+  void reset(const casa::uInt size=0, const casa::uInt tsysSize=0);
+  void initialize(const casa::uInt size, const casa::uInt tsysSize);
   /**
     * check the initialization state 
     */ 
   casa::Bool state() const;
+  /**
+    * replace NaN values with (normal) values at the same channels in the given spetrum.
+    * (CAS-2776; 2011/04/07 by Wataru Kawasaki)
+    */
+  void replaceNaN();
 
 private:
-  void addSpectrum( const casa::Vector<casa::Float>& v,
-                    const casa::Vector<casa::Bool>& m,
-                    casa::Float weight);
-
-  casa::Float addTsys(const casa::Vector<casa::Float>& v);
-  casa::Float addInterval(casa::Double inter);
-  void addTime(casa::Double t);
+  void addSpectrum(const casa::Vector<casa::Float>& v,
+                   const casa::Vector<casa::Bool>& m,
+                   const casa::Vector<casa::Float>& tsys,
+                   const casa::Double interval,
+                   const casa::Double time);
+  void doAddSpectrum(const casa::Vector<casa::Float>& v,
+                     const casa::Vector<casa::Bool>& m,
+                     const casa::Vector<casa::Float>& tsys,
+                     const casa::Double interval,
+                     const casa::Double time,
+                     const casa::Bool ignoreMask);
+  casa::Float getTotalWeight(const casa::MaskedArray<casa::Float>& data,
+                             const casa::Vector<casa::Float>& tsys,
+                             const casa::Double interval,
+                             const casa::Double time,
+                             const casa::Bool ignoreMask);
+  casa::Float addTsys(const casa::Vector<casa::Float>& v, casa::Bool ignoreMask);
+  casa::Float addInterval(casa::Double inter, casa::Bool ignoreMask);
+  void addTime(casa::Double t, casa::Bool ignoreMask);
 
   WeightType weightType_;
@@ -106 +125 @@
   casa::MaskedArray<casa::uInt> n_;
 
+  //these three are used for normalise() (CAS-2776; 2011/04/07 by WK)
+  casa::MaskedArray<casa::Float> spectrumNoMask_;
+  casa::MaskedArray<casa::Float> weightSumNoMask_;
+  casa::MaskedArray<casa::uInt> nNoMask_;
+
   casa::Vector<casa::Bool> userMask_;
 
-  casa::Vector<casa::Float> tsysSum_;
-  casa::Double timeSum_;
-  casa::Double intervalSum_;
+  casa::Vector<casa::Float> tsysSum_, tsysSumNoMask_;
+  casa::Double timeSum_, timeSumNoMask_;
+  casa::Double intervalSum_, intervalSumNoMask_;
 };
 

branches/casa-prerelease/pre-asap/src/STMath.cpp

@@ -207 +207 @@
       ROScalarColumn<Double> tmp(tin, "TIME");
       Double td;tmp.get(0,td);
-      Table basesubt = tin(tin.col("BEAMNO") == Int(rec.asuInt("BEAMNO"))
-                       && tin.col("IFNO") == Int(rec.asuInt("IFNO"))
-                       && tin.col("POLNO") == Int(rec.asuInt("POLNO")) );
+      Table basesubt = tin( tin.col("BEAMNO") == Int(rec.asuInt("BEAMNO"))
+                         && tin.col("IFNO") == Int(rec.asuInt("IFNO"))
+                         && tin.col("POLNO") == Int(rec.asuInt("POLNO")) );
       Table subt;
       if ( avmode == "SOURCE") {
-        subt = basesubt( basesubt.col("SRCNAME") == rec.asString("SRCNAME") );
+        subt = basesubt( basesubt.col("SRCNAME") == rec.asString("SRCNAME"));
       } else if (avmode == "SCAN") {
-        //subt = basesubt( basesubt.col("SCANNO") == Int(rec.asuInt("SCANNO")) );
-        subt = basesubt( basesubt.col("SCANNO") == Int(rec.asuInt("SCANNO"))
-                         && basesubt.col("SRCNAME") == rec.asString("SRCNAME") );
+        subt = basesubt( basesubt.col("SRCNAME") == rec.asString("SRCNAME") 
+                      && basesubt.col("SCANNO") == Int(rec.asuInt("SCANNO")) );
       } else {
         subt = basesubt;
@@ -266 +265 @@
         acc.add(spec, !bflag, tsys, inter, time);
       }
+
+      // If there exists a channel at which all the input spectra are masked, 
+      // spec has 'nan' values for that channel and it may affect the following 
+      // processes. To avoid this, replacing 'nan' values in spec with 
+      // weighted-mean of all spectra in the following line. 
+      // (done for CAS-2776, 2011/04/07 by Wataru Kawasaki)
+      acc.replaceNaN();
     }
     const Vector<Bool>& msk = acc.getMask();
@@ -278 +284 @@
       Vector<uChar> flg(msk.shape());
       convertArray(flg, !msk);
+      for (uInt k = 0; k < flg.nelements(); ++k) {
+        uChar userFlag = 1 << 7;
+        if (msk[k]==True) userFlag = 0 << 7;
+        flg(k) = userFlag;
+      }
+
       flagColOut.put(i, flg);
       specColOut.put(i, acc.getSpectrum());
@@ -627 +639 @@
   // check if vector size is equal to nrow
   Vector<Float> fact( val ) ;
-  if ( fact.nelements() != in->nrow() ) {
+  if (fact.nelements() != uInt(in->nrow())) {
     throw( AipsError("Vector size must be 1 or be same as number of row.") ) ;
   }
@@ -688 +700 @@
   // some checks 
   vector<uInt> nchans;
-  for ( uInt i = 0 ; i < in->nrow() ; i++ ) {
-    nchans.push_back( (in->getSpectrum( i )).size() ) ;
+  for (Int i = 0 ; i < in->nrow() ; i++) {
+    nchans.push_back((in->getSpectrum(i)).size());
   }
   //Vector<uInt> mchans( nchans ) ;
@@ -2248 +2260 @@
   Table outtab = out->table();
 
-  const uInt ntau = uInt(tau.size());
+  const Int ntau = uInt(tau.size());
   std::vector<float>::const_iterator tauit = tau.begin();
   AlwaysAssert((ntau == 1 || ntau == in->nif() || ntau == in->nif() * in->npol()),
@@ -2292 +2304 @@
 {
   CountedPtr< Scantable > out = getScantable(in, false);
-  Table& table = out->table();
-  ArrayColumn<Float> specCol(table, "SPECTRA");
-  ArrayColumn<uChar> flagCol(table, "FLAGTRA");
-  Vector<Float> spec;
-  Vector<uChar> flag;
-  for ( uInt i=0; i<table.nrow(); ++i) {
-    specCol.get(i, spec);
-    flagCol.get(i, flag);
-    Vector<Bool> mask(flag.nelements());
-    convertArray(mask, flag);
-    Vector<Float> specout;
-    Vector<Bool> maskout;
-    if ( kernel == "hanning" ) {
-      mathutil::hanning(specout, maskout, spec , !mask);
-      convertArray(flag, !maskout);
-    } else if (  kernel == "rmedian" ) {
-      mathutil::runningMedian(specout, maskout, spec , mask, width);
-      convertArray(flag, maskout);
-    } else if ( kernel == "poly" ) {
-      mathutil::polyfit(specout, maskout, spec, !mask, width, order);
-      convertArray(flag, !maskout);
-    }
-    flagCol.put(i, flag);
-    specCol.put(i, specout);
+  Table table = out->table();
+
+  TableIterator iter(table, "IFNO");
+  while (!iter.pastEnd()) {
+    Table tab = iter.table();
+    ArrayColumn<Float> specCol(tab, "SPECTRA");
+    ArrayColumn<uChar> flagCol(tab, "FLAGTRA");
+    Vector<Float> spec;
+    Vector<uChar> flag;
+    for (uInt i = 0; i < tab.nrow(); ++i) {
+      specCol.get(i, spec);
+      flagCol.get(i, flag);
+      Vector<Bool> mask(flag.nelements());
+      convertArray(mask, flag);
+      Vector<Float> specout;
+      Vector<Bool> maskout;
+      if (kernel == "hanning") {
+        mathutil::hanning(specout, maskout, spec, !mask);
+        convertArray(flag, !maskout);
+      } else if (kernel == "rmedian") {
+        mathutil::runningMedian(specout, maskout, spec , mask, width);
+        convertArray(flag, maskout);
+      } else if (kernel == "poly") {
+        mathutil::polyfit(specout, maskout, spec, !mask, width, order);
+        convertArray(flag, !maskout);
+      }
+
+      for (uInt j = 0; j < flag.nelements(); ++j) {
+        uChar userFlag = 1 << 7;
+        if (maskout[j]==True) userFlag = 0 << 7;
+        flag(j) = userFlag;
+      }
+
+      flagCol.put(i, flag);
+      specCol.put(i, specout);
+    }
+  ++iter;
   }
   return out;
@@ -2822 +2847 @@
         if (fstart-1 > 0 ) {
           interp = spec[fstart-1];
-          if (fend+1 < spec.nelements()) {
+          if (fend+1 < Int(spec.nelements())) {
             interp = (interp+spec[fend+1])/2.0;
           }

branches/casa-prerelease/pre-asap/src/python_Scantable.cpp

@@ -155 +155 @@
     .def("_getflagtrafast", &ScantableWrapper::getFlagtraFast,
          (boost::python::arg("whichrow")=0) )
-    //.def("_sspline_baseline", &ScantableWrapper::smoothingSplineBaseline,
-    // (boost::python::arg("thres")=3.0,
-    //  boost::python::arg("niter")=1) )
+    //.def("_sspline_baseline", &ScantableWrapper::smoothingSplineBaseline)
     //.def("_test_cin", &ScantableWrapper::testCin)
   ;