Changeset 779 for branches/Release12

Timestamp:

12/07/05 14:14:45 (19 years ago)

Author:

mar637

Message:

log message issues in calcazel

Location:

branches/Release12/src

Files:

• SDMath.cc (modified) (89 diffs)
• SDMemTable.cc (modified) (4 diffs)

branches/Release12/src/SDMath.cc

@@ -101 +101 @@
 }
 
-SDMath& SDMath::operator=(const SDMath& other) 
+SDMath& SDMath::operator=(const SDMath& other)
 {
   if (this != &other) {
@@ -113 +113 @@
 
 
-SDMemTable* SDMath::frequencyAlignment(const SDMemTable& in, 
-                                       const String& refTime, 
-                                       const String& method, 
-                                       Bool perFreqID)
+SDMemTable* SDMath::frequencyAlignment(const SDMemTable& in,
+                                       const String& refTime,
+                                       const String& method,
+                                       Bool perFreqID)
 {
   // Get frame info from Table
@@ -136 +136 @@
 
 
-CountedPtr<SDMemTable> 
+CountedPtr<SDMemTable>
 SDMath::average(const std::vector<CountedPtr<SDMemTable> >& in,
-                const Vector<Bool>& mask, Bool scanAv,
-                const String& weightStr, Bool alignFreq)
+                const Vector<Bool>& mask, Bool scanAv,
+                const String& weightStr, Bool alignFreq)
 // Weighted averaging of spectra from one or more Tables.
 {
-  // Convert weight type  
+  // Convert weight type
   WeightType wtType = NONE;
   convertWeightString(wtType, weightStr, True);
@@ -242 +242 @@
     // Should check that the frequency tables don't change if doing
     // FreqAlignment
-    
+
     // Attach columns to Table
      const Table& tabIn = in[iTab]->table();
@@ -258 +258 @@
         IPosition shp2 = in[iTab]->rowAsMaskedArray(iRow).shape();
         if (!shp.isEqual(shp2)) {
-          delete pTabOut;
+          delete pTabOut;
            throw (AipsError("Shapes for all rows must be the same"));
         }
 
-        // If we are not doing scan averages, make checks for source
-        // and frequency setup and warn if averaging across them
+        // If we are not doing scan averages, make checks for source
+        // and frequency setup and warn if averaging across them
         scanIDCol.getScalar(iRow, scanID);
 
-        // Get quantities from columns
+        // Get quantities from columns
         srcNameCol.getScalar(iRow, sourceName);
         mjdCol.get(iRow, time);
@@ -273 +273 @@
         fqIDCol.get(iRow, freqID);
 
-        // Initialize first source and freqID
+        // Initialize first source and freqID
         if (iRow==0 && iTab==0) {
           sourceNameStart = sourceName;
@@ -279 +279 @@
         }
 
-        // If we are doing scan averages, see if we are at the end of
-        // an accumulation period (scan).  We must check soutce names
-        // too, since we might have two tables with one scan each but
-        // different source names; we shouldn't average different
-        // sources together
+        // If we are doing scan averages, see if we are at the end of
+        // an accumulation period (scan).  We must check soutce names
+        // too, since we might have two tables with one scan each but
+        // different source names; we shouldn't average different
+        // sources together
         if (scanAv && ( (scanID != oldScanID)  ||
                         (iRow==0 && iTab>0 && sourceName!=oldSourceName))) {
 
-          // Normalize data in 'sum' accumulation array according to
-          // weighting scheme
-           normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType, 
-                     asap::ChanAxis, nAxesSub);
-
-           // Get ScanContainer for the first row of this averaged Scan
+          // Normalize data in 'sum' accumulation array according to
+          // weighting scheme
+           normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType,
+                     asap::ChanAxis, nAxesSub);
+
+           // Get ScanContainer for the first row of this averaged Scan
            SDContainer scOut = in[iTab]->getSDContainer(rowStart);
 
-           // Fill scan container. The source and freqID come from the
-           // first row of the first table that went into this average
-           // ( should be the same for all rows in the scan average)
+           // Fill scan container. The source and freqID come from the
+           // first row of the first table that went into this average
+           // ( should be the same for all rows in the scan average)
 
            Float nR(nAccum);
            fillSDC(scOut, sum.getMask(), sum.getArray(), tSysSum/nR, outScanID,
-                   timeSum/nR, intSum, sourceNameStart, freqIDStart);
-           
-           // Write container out to Table
+                   timeSum/nR, intSum, sourceNameStart, freqIDStart);
+
+           // Write container out to Table
            pTabOut->putSDContainer(scOut);
-           
-           // Reset accumulators           
+
+           // Reset accumulators
            sum = 0.0;
            sumSq = 0.0;
@@ -315 +315 @@
            timeSum = 0.0;
            intSum = 0.0;
-           nPts = 0.0;
-
-           // Increment
+           nPts = 0.0;
+
+           // Increment
            rowStart = iRow;              // First row for next accumulation
            tableStart = iTab;            // First table for next accumulation
            sourceNameStart = sourceName; // First source name for next
-                                         // accumulation
+                                         // accumulation
            freqIDStart = freqID;         // First FreqID for next accumulation
 
            oldScanID = scanID;
            outScanID += 1;               // Scan ID for next
-                                         // accumulation period
+                                         // accumulation period
         }
 
-        // Accumulate
-        accumulate(timeSum, intSum, nAccum, sum, sumSq, nPts, 
-                   tSysSum, tSysSqSum, tSys, 
-                   nInc, mask, time, interval, in, iTab, iRow, asap::ChanAxis, 
+        // Accumulate
+        accumulate(timeSum, intSum, nAccum, sum, sumSq, nPts,
+                   tSysSum, tSysSqSum, tSys,
+                   nInc, mask, time, interval, in, iTab, iRow, asap::ChanAxis,
                    nAxesSub, useMask, wtType);
-        oldSourceName = sourceName;
-        oldFreqID = freqID;
+        oldSourceName = sourceName;
+        oldFreqID = freqID;
      }
   }
@@ -347 +347 @@
   // scheme
 
-  normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType, 
-            asap::ChanAxis, nAxesSub);
+  normalize(sum, sumSq, tSysSqSum, nPts, intSum, wtType,
+            asap::ChanAxis, nAxesSub);
 
   // Create and fill container.  The container we clone will be from
@@ -356 +356 @@
   SDContainer scOut = in[tableStart]->getSDContainer(rowStart);
   fillSDC(scOut, sum.getMask(), sum.getArray(), tSysSum/nR, outScanID,
-          timeSum/nR, intSum, sourceNameStart, freqIDStart);
+          timeSum/nR, intSum, sourceNameStart, freqIDStart);
   pTabOut->putSDContainer(scOut);
   pTabOut->resetCursor();
@@ -365 +365 @@
 
 
-CountedPtr<SDMemTable> 
+CountedPtr<SDMemTable>
 SDMath::binaryOperate(const CountedPtr<SDMemTable>& left,
-                      const CountedPtr<SDMemTable>& right,
-                      const String& op, Bool preserve, Bool doTSys)
+                      const CountedPtr<SDMemTable>& right,
+                      const String& op, Bool preserve, Bool doTSys)
 {
 
@@ -399 +399 @@
   }
 
-  // Input Tables 
+  // Input Tables
   const Table& tLeft = left->table();
   const Table& tRight = right->table();
@@ -413 +413 @@
   Array<Float> tSysLeftArr, tSysRightArr;
   if (doTSys) tSysRightCol.get(0, tSysRightArr);
-  MaskedArray<Float>* pMRight = 
+  MaskedArray<Float>* pMRight =
     new MaskedArray<Float>(right->rowAsMaskedArray(0));
 
@@ -424 +424 @@
   // Loop over rows
   for (uInt i=0; i<nRowLeft; i++) {
-    
-    // Get data 
+
+    // Get data
     MaskedArray<Float> mLeft(left->rowAsMaskedArray(i));
     IPosition shpLeft = mLeft.shape();
     if (doTSys) tSysLeftCol.get(i, tSysLeftArr);
-    
+
     if (nRowRight>1) {
       delete pMRight;
@@ -444 +444 @@
     if (doTSys) {
       if (!tSysRightArr.shape().isEqual(tSysRightArr.shape())) {
-        delete pTabOut;
-        delete pMRight;
-        throw(AipsError("left and right Tsys data are not conformant"));
+        delete pTabOut;
+        delete pMRight;
+        throw(AipsError("left and right Tsys data are not conformant"));
       }
       if (!shpRight.isEqual(tSysRightArr.shape())) {
-        delete pTabOut;
-        delete pMRight;
-        throw(AipsError("left and right scan tables are not conformant"));
+        delete pTabOut;
+        delete pMRight;
+        throw(AipsError("left and right scan tables are not conformant"));
       }
     }
@@ -459 +459 @@
 
      // Operate on data and TSys
-     if (what==0) {                               
+     if (what==0) {
         MaskedArray<Float> tmp = mLeft + *pMRight;
         putDataInSDC(sc, tmp.getArray(), tmp.getMask());
@@ -476 +476 @@
         if (doTSys) sc.putTsys(tSysLeftArr/tSysRightArr);
      } else if (what==4) {
-       if (preserve) {     
-         MaskedArray<Float> tmp = (tSysRightArr * mLeft / *pMRight) - 
-           tSysRightArr;
-         putDataInSDC(sc, tmp.getArray(), tmp.getMask());
+       if (preserve) {
+         MaskedArray<Float> tmp = (tSysRightArr * mLeft / *pMRight) -
+           tSysRightArr;
+         putDataInSDC(sc, tmp.getArray(), tmp.getMask());
        } else {
-         MaskedArray<Float> tmp = (tSysRightArr * mLeft / *pMRight) - 
-           tSysLeftArr;
-         putDataInSDC(sc, tmp.getArray(), tmp.getMask());
+         MaskedArray<Float> tmp = (tSysRightArr * mLeft / *pMRight) -
+           tSysLeftArr;
+         putDataInSDC(sc, tmp.getArray(), tmp.getMask());
        }
        sc.putTsys(tSysRightArr);
@@ -493 +493 @@
   if (pMRight) delete pMRight;
   pTabOut->resetCursor();
-  
+
   return CountedPtr<SDMemTable>(pTabOut);
 }
@@ -499 +499 @@
 
 std::vector<float> SDMath::statistic(const CountedPtr<SDMemTable>& in,
-                                     const Vector<Bool>& mask,
-                                     const String& which, Int row) const
+                                     const Vector<Bool>& mask,
+                                     const String& which, Int row) const
 //
 // Perhaps iteration over pol/beam/if should be in here
@@ -519 +519 @@
   uInt iStart = 0;
   uInt iEnd = in->nRow()-1;
-//  
+//
   if (row>=0) {
      iStart = row;
@@ -569 +569 @@
     CoordinateSystem cSys;
     cSys.addCoordinate(in.getSpectralCoordinate(j));
-    CoordinateSystem cSysBin = 
+    CoordinateSystem cSysBin =
       CoordinateUtil::makeBinnedCoordinateSystem(factors, cSys, False);
 //
@@ -584 +584 @@
 
 // Loop over rows and bin along channel axis
-  
+
   for (uInt i=0; i < in.nRow(); ++i) {
     SDContainer sc = in.getSDContainer(i);
@@ -603 +603 @@
     putDataInSDC(sc, marrout.getArray(), marrout.getMask());
 
-// Bin up Tsys.  
+// Bin up Tsys.
 
     Array<Bool> allGood(tSys.shape(),True);
     MaskedArray<Float> tSysIn(tSys, allGood, True);
 //
-    MaskedArray<Float> tSysOut;    
+    MaskedArray<Float> tSysOut;
     LatticeUtilities::bin(tSysOut, tSysIn, asap::ChanAxis, width);
     sc.putTsys(tSysOut.getArray());
@@ -618 +618 @@
 
 SDMemTable* SDMath::resample(const SDMemTable& in, const String& methodStr,
-                             Float width)
+                             Float width)
 //
 // Should add the possibility of width being specified in km/s. This means
-// that for each freqID (SpectralCoordinate) we will need to convert to an 
-// average channel width (say at the reference pixel).  Then we would need  
-// to be careful to make sure each spectrum (of different freqID) 
+// that for each freqID (SpectralCoordinate) we will need to convert to an
+// average channel width (say at the reference pixel).  Then we would need
+// to be careful to make sure each spectrum (of different freqID)
 // is the same length.
 //
@@ -632 +632 @@
          SpectralCoordinate sC = in.getSpectralCoordinate(j);
       }
-   } 
+   }
 
 // Interpolation method
@@ -691 +691 @@
 
   Array<Float> valuesOut;
-  Array<Bool> maskOut;  
+  Array<Bool> maskOut;
   Array<Float> tSysOut;
   Array<Bool> tSysMaskIn(shapeIn,True);
@@ -702 +702 @@
 
 // Get data and Tsys
-    
+
      const Array<Float>& tSysIn = sc.getTsys();
      const MaskedArray<Float>& dataIn(in.rowAsMaskedArray(i));
@@ -710 +710 @@
 // Interpolate data
 
-     InterpolateArray1D<Float,Float>::interpolate(valuesOut, maskOut, xOut, 
+     InterpolateArray1D<Float,Float>::interpolate(valuesOut, maskOut, xOut,
                                                   xIn, valuesIn, maskIn,
                                                   interpMethod, True, True);
@@ -718 +718 @@
 // Interpolate TSys
 
-     InterpolateArray1D<Float,Float>::interpolate(tSysOut, tSysMaskOut, xOut, 
+     InterpolateArray1D<Float,Float>::interpolate(tSysOut, tSysMaskOut, xOut,
                                                   xIn, tSysIn, tSysMaskIn,
                                                   interpMethod, True, True);
@@ -739 +739 @@
    SDMemTable* pOut = new SDMemTable(in,False);
    const Table& tOut = pOut->table();
-   ArrayColumn<Float> specCol(tOut,"SPECTRA");  
-   ArrayColumn<Float> tSysCol(tOut,"TSYS");  
+   ArrayColumn<Float> specCol(tOut,"SPECTRA");
+   ArrayColumn<Float> tSysCol(tOut,"TSYS");
    Array<Float> tSysArr;
 
@@ -797 +797 @@
   const Table& tabIn = in.table();
 
-// Shape of input and output data 
+// Shape of input and output data
 
   const IPosition& shapeIn = in.rowAsMaskedArray(0).shape();
@@ -826 +826 @@
      tSysCol.attach(tabIn,"TSYS");
   }
-// 
+//
   const Bool useMask = (mask.nelements() == shapeIn(asap::ChanAxis));
 
@@ -838 +838 @@
       Array<Float>& arr = marr.getRWArray();
       const Array<Bool>& barr = marr.getMask();
-      
+
 // Get Tsys
 
@@ -852 +852 @@
       ReadOnlyArrayIterator<Bool> itMaskPlane(barr, axes);
       ReadOnlyArrayIterator<Float>* pItTsysPlane = 0;
-      if (wtType==TSYS) 
-        pItTsysPlane = new ReadOnlyArrayIterator<Float>(tSys, axes);
+      if (wtType==TSYS)
+        pItTsysPlane = new ReadOnlyArrayIterator<Float>(tSys, axes);
 
 // Accumulations
@@ -868 +868 @@
         Vector<Float> t1(nChan); t1 = 0.0;
         Vector<Bool> t2(nChan); t2 = True;
-        Float tSys = 0.0;
+        Float tSys = 0.0;
         MaskedArray<Float> vecSum(t1,t2);
         Float norm = 0.0;
@@ -876 +876 @@
 //
            ReadOnlyVectorIterator<Float>* pItTsysVec = 0;
-           if (wtType==TSYS) {
-              pItTsysVec = 
-                new ReadOnlyVectorIterator<Float>(pItTsysPlane->array(), 1);
-           }             
-//
-           while (!itDataVec.pastEnd()) {     
+           if (wtType==TSYS) {
+              pItTsysVec =
+                new ReadOnlyVectorIterator<Float>(pItTsysPlane->array(), 1);
+           }
+//
+           while (!itDataVec.pastEnd()) {
 
 // Create MA of data & mask (optionally including OTF mask) and  get variance for this spectrum
@@ -887 +887 @@
               if (useMask) {
                  const MaskedArray<Float> spec(itDataVec.vector(),
-                                               mask&&itMaskVec.vector());
+                                               mask&&itMaskVec.vector());
                  if (wtType==VAR) {
                     fac = 1.0 / variance(spec);
                  } else if (wtType==TSYS) {
-                    tSys = pItTsysVec->vector()[0];      // Drop pseudo channel dependency
+                    tSys = pItTsysVec->vector()[0];      // Drop pseudo channel dependency
                     fac = 1.0 / tSys / tSys;
-                 }                    
+                 }
               } else {
                  const MaskedArray<Float> spec(itDataVec.vector(),
-                                               itMaskVec.vector());
+                                               itMaskVec.vector());
                  if (wtType==VAR) {
                     fac = 1.0 / variance(spec);
                  } else if (wtType==TSYS) {
-                    tSys = pItTsysVec->vector()[0];      // Drop pseudo channel dependency
+                    tSys = pItTsysVec->vector()[0];      // Drop pseudo channel dependency
                     fac = 1.0 / tSys / tSys;
                  }
@@ -908 +908 @@
 
               const MaskedArray<Float> spec(fac*itDataVec.vector(),
-                                            itMaskVec.vector());
+                                            itMaskVec.vector());
               vecSum += spec;
               norm += fac;
@@ -916 +916 @@
               itDataVec.next();
               itMaskVec.next();
-              if (wtType==TSYS) pItTsysVec->next();
+              if (wtType==TSYS) pItTsysVec->next();
            }
-           
+
 // Clean up
 
-           if (pItTsysVec) {
-              delete pItTsysVec; 
+           if (pItTsysVec) {
+              delete pItTsysVec;
               pItTsysVec = 0;
-           }            
+           }
         }
 
@@ -936 +936 @@
         IPosition pos = itDataPlane.pos();
 
-// Write out data. This is a bit messy. We have to reform the Vector 
+// Write out data. This is a bit messy. We have to reform the Vector
 // accumulator into an Array of shape (1,1,1,nChan)
 
         start = pos;
-        end = pos; 
+        end = pos;
         end(asap::ChanAxis) = nChan-1;
         outData(start,end) = vecSum.getArray().reform(vecShapeOut);
@@ -949 +949 @@
         itDataPlane.next();
         itMaskPlane.next();
-        if (wtType==TSYS) pItTsysPlane->next();
+        if (wtType==TSYS) pItTsysPlane->next();
       }
 
@@ -974 +974 @@
 
 
-SDMemTable* SDMath::smooth(const SDMemTable& in, 
-                           const casa::String& kernelType,
-                           casa::Float width, Bool doAll)
+SDMemTable* SDMath::smooth(const SDMemTable& in,
+                           const casa::String& kernelType,
+                           casa::Float width, Bool doAll)
 //
 // Should smooth TSys as well
@@ -1014 +1014 @@
       Array<Bool> maskIn = dataIn.getMask();
 
-      Array<Float> valuesIn2 = valuesIn(start,end);       // ref to valuesIn 
+      Array<Float> valuesIn2 = valuesIn(start,end);       // ref to valuesIn
       Array<Bool> maskIn2 = maskIn(start,end);
 
@@ -1021 +1021 @@
       VectorIterator<Bool> itMask(maskIn2, asap::ChanAxis);
       while (!itValues.pastEnd()) {
-        
-        // Smooth
-        if (kernelType==VectorKernel::HANNING) {
-          mathutil::hanning(valuesOut, maskOut, itValues.vector(), 
-                            itMask.vector());
-          itMask.vector() = maskOut;
-        } else {
-          mathutil::replaceMaskByZero(itValues.vector(), itMask.vector());
-          conv.linearConv(valuesOut, itValues.vector());
-        }
-
-        itValues.vector() = valuesOut;
-        itValues.next();
-        itMask.next();
+
+        // Smooth
+        if (kernelType==VectorKernel::HANNING) {
+          mathutil::hanning(valuesOut, maskOut, itValues.vector(),
+                            itMask.vector());
+          itMask.vector() = maskOut;
+        } else {
+          mathutil::replaceMaskByZero(itValues.vector(), itMask.vector());
+          conv.linearConv(valuesOut, itValues.vector());
+        }
+
+        itValues.vector() = valuesOut;
+        itValues.next();
+        itMask.next();
       }
-      
+
       // Create and put back
       SDContainer sc = in.getSDContainer(ri);
@@ -1049 +1049 @@
 
 
-SDMemTable* SDMath::convertFlux(const SDMemTable& in, Float D, Float etaAp, 
-                                Float JyPerK, Bool doAll)
-// 
+SDMemTable* SDMath::convertFlux(const SDMemTable& in, Float D, Float etaAp,
+                                Float JyPerK, Bool doAll)
+//
 // etaAp = aperture efficiency (-1 means find)
 // D     = geometric diameter (m)  (-1 means find)
-// JyPerK 
+// JyPerK
 //
 {
@@ -1063 +1063 @@
   // mJy or mK) Also automatically find out what we are converting to
   // according to the flux unit
-  Unit fluxUnit(sh.fluxunit); 
+  Unit fluxUnit(sh.fluxunit);
   Unit K(String("K"));
   Unit JY(String("Jy"));
 
   Bool toKelvin = True;
-  Double cFac = 1.0;    
+  Double cFac = 1.0;
 
   if (fluxUnit==JY) {
@@ -1075 +1075 @@
     Quantum<Double> t2 = t.get(JY);
     cFac = (t2 / t).getValue();               // value to Jy
-    
+
     toKelvin = True;
     sh.fluxunit = "K";
@@ -1083 +1083 @@
     Quantum<Double> t2 = t.get(K);
     cFac = (t2 / t).getValue();              // value to K
-    
+
     toKelvin = False;
     sh.fluxunit = "Jy";
@@ -1091 +1091 @@
 
   pTabOut->putSDHeader(sh);
-  
+
   // Make sure input values are converted to either Jy or K first...
   Float factor = cFac;
@@ -1121 +1121 @@
     scaleByVector(pTabOut, in, doAll, factors, False);
   } else {
-    
+
     // OK now we must deal with automatic look up of values.
     // We must also deal with the fact that the factors need
     // to be computed per IF and may be different and may
     // change per integration.
-    
+
     pushLog("Looking up conversion factors");
     convertBrightnessUnits (pTabOut, in, toKelvin, cFac, doAll);
@@ -1134 +1134 @@
 
 
-SDMemTable* SDMath::gainElevation(const SDMemTable& in, 
-                                  const Vector<Float>& coeffs,
-                                  const String& fileName,
-                                  const String& methodStr, Bool doAll)
+SDMemTable* SDMath::gainElevation(const SDMemTable& in,
+                                  const Vector<Float>& coeffs,
+                                  const String& fileName,
+                                  const String& methodStr, Bool doAll)
 {
 
@@ -1149 +1149 @@
   Vector<Float> x = elev.getColumn();
   x *= Float(180 / C::pi);                        // Degrees
-  
+
   const uInt nC = coeffs.nelements();
   if (fileName.length()>0 && nC>0) {
     throw(AipsError("You must choose either polynomial coefficients or an ascii file, not both"));
   }
-  
+
   // Correct
   if (nC>0 || fileName.length()==0) {
@@ -1160 +1160 @@
      Bool throwIt = True;
      Instrument inst = SDAttr::convertInstrument (sh.antennaname, throwIt);
-     
+
      // Set polynomial
      Polynomial<Float>* pPoly = 0;
@@ -1175 +1175 @@
        msg = String("built in");
      }
-     
+
      if (coeff.nelements()>0) {
        pPoly->setCoefficients(coeff);
@@ -1182 +1182 @@
        throw(AipsError("There is no known gain-elevation polynomial known for this instrument"));
      }
-     pushLog("Making polynomial correction with "+msg+" coefficients");
+     ostringstream oss;
+     oss << "Making polynomial correction with " << msg << " coefficients:" << endl;
+     oss << "   " <<  coeff;
+     pushLog(String(oss));
      const uInt nRow = in.nRow();
      Vector<Float> factor(nRow);
@@ -1192 +1195 @@
 
   } else {
-    
+
     // Indicate which columns to read from ascii file
     String col0("ELEVATION");
     String col1("FACTOR");
-    
+
     // Read and correct
-    
+
     pushLog("Making correction from ascii Table");
-    scaleFromAsciiTable (pTabOut, in, fileName, col0, col1, 
-                         methodStr, doAll, x, True);
+    scaleFromAsciiTable (pTabOut, in, fileName, col0, col1,
+                         methodStr, doAll, x, True);
   }
 
   return pTabOut;
 }
- 
+
 
 SDMemTable* SDMath::opacity(const SDMemTable& in, Float tau, Bool doAll)
@@ -1256 +1259 @@
       throw(AipsError("Only linear polarizations are supported"));
    }
-//    
+//
    const Table& tabIn = in.table();
-   ArrayColumn<Float> specCol(tabIn,"SPECTRA");  
+   ArrayColumn<Float> specCol(tabIn,"SPECTRA");
    IPosition start(asap::nAxes,0);
    IPosition end(asap::nAxes);
 
 // Set cursor slice. Assumes shape the same for all rows
- 
+
    setCursorSlice (start, end, doAll, in);
    IPosition start3(start);
    start3(asap::PolAxis) = 2;                 // Real(XY)
    IPosition end3(end);
-   end3(asap::PolAxis) = 2;   
+   end3(asap::PolAxis) = 2;
 //
    IPosition start4(start);
@@ -1274 +1277 @@
    IPosition end4(end);
    end4(asap::PolAxis) = 3;
-//  
+//
    uInt nRow = in.nRow();
    Array<Float> data;
@@ -1280 +1283 @@
       specCol.get(i,data);
       IPosition shape = data.shape();
- 
+
       // Get polarization slice references
       Array<Float> C3 = data(start3,end3);
       Array<Float> C4 = data(start4,end4);
-    
+
       // Rotate
       SDPolUtil::rotatePhase(C3, C4, value);
-   
+
       // Put
       specCol.put(i,data);
    }
-}     
+}
 
 
@@ -1310 +1313 @@
       throw(AipsError("Only linear polarizations are supported"));
    }
-//    
+//
    const Table& tabIn = in.table();
-   ArrayColumn<Float> specCol(tabIn,"SPECTRA");  
-   ROArrayColumn<Float> stokesCol(tabIn,"STOKES");  
+   ArrayColumn<Float> specCol(tabIn,"SPECTRA");
+   ROArrayColumn<Float> stokesCol(tabIn,"STOKES");
    IPosition start(asap::nAxes,0);
    IPosition end(asap::nAxes);
 
 // Set cursor slice. Assumes shape the same for all rows
- 
+
    setCursorSlice (start, end, doAll, in);
 //
@@ -1324 +1327 @@
    start1(asap::PolAxis) = 0;                // C1 (XX)
    IPosition end1(end);
-   end1(asap::PolAxis) = 0;   
+   end1(asap::PolAxis) = 0;
 //
    IPosition start2(start);
    start2(asap::PolAxis) = 1;                 // C2 (YY)
    IPosition end2(end);
-   end2(asap::PolAxis) = 1;   
+   end2(asap::PolAxis) = 1;
 //
    IPosition start3(start);
    start3(asap::PolAxis) = 2;                 // C3 ( Real(XY) )
    IPosition end3(end);
-   end3(asap::PolAxis) = 2;   
+   end3(asap::PolAxis) = 2;
 //
    IPosition startI(start);
    startI(asap::PolAxis) = 0;                 // I
    IPosition endI(end);
-   endI(asap::PolAxis) = 0;   
+   endI(asap::PolAxis) = 0;
 //
    IPosition startQ(start);
    startQ(asap::PolAxis) = 1;                 // Q
    IPosition endQ(end);
-   endQ(asap::PolAxis) = 1;   
+   endQ(asap::PolAxis) = 1;
 //
    IPosition startU(start);
    startU(asap::PolAxis) = 2;                 // U
    IPosition endU(end);
-   endU(asap::PolAxis) = 2;   
+   endU(asap::PolAxis) = 2;
 
 //
@@ -1358 +1361 @@
       stokesCol.get(i,stokes);
       IPosition shape = data.shape();
- 
+
 // Get linear polarization slice references
-  
+
       Array<Float> C1 = data(start1,end1);
       Array<Float> C2 = data(start2,end2);
@@ -1370 +1373 @@
       Array<Float> Q = stokes(startQ,endQ);
       Array<Float> U = stokes(startU,endU);
-    
+
 // Rotate
- 
+
       SDPolUtil::rotateLinPolPhase(C1, C2, C3, I, Q, U, value);
-   
+
 // Put
-    
+
       specCol.put(i,data);
    }
-}     
+}
 
 // 'private' functions
 
-void SDMath::convertBrightnessUnits (SDMemTable* pTabOut, const SDMemTable& in, 
+void SDMath::convertBrightnessUnits (SDMemTable* pTabOut, const SDMemTable& in,
                                      Bool toKelvin, Float cFac, Bool doAll)
 {
@@ -1438 +1441 @@
 // Loop over rows and apply correction factor
 
-   Float factor = 1.0;  
+   Float factor = 1.0;
    const uInt axis = asap::ChanAxis;
    for (uInt i=0; i < in.nRow(); ++i) {
@@ -1457 +1460 @@
 
 // Now the conversion factor depends only upon frequency
-// So we need to iterate through by IF only giving 
+// So we need to iterate through by IF only giving
 // us BEAM/POL/CHAN cubes
 
@@ -1479 +1482 @@
 SDMemTable* SDMath::frequencyAlign(const SDMemTable& in,
                                    MFrequency::Types freqSystem,
-                                   const String& refTime, 
+                                   const String& refTime,
                                    const String& methodStr,
                                    Bool perFreqID)
@@ -1501 +1504 @@
    Vector<Double> times = mjdCol.getColumn();
 
-// Generate DataDesc table 
- 
+// Generate DataDesc table
+
    Matrix<uInt> ddIdx;
    SDDataDesc dDesc;
-   generateDataDescTable(ddIdx, dDesc, nIF, in, tabIn, srcCol, 
-                          fqIDCol, perFreqID);
+   generateDataDescTable(ddIdx, dDesc, nIF, in, tabIn, srcCol,
+                          fqIDCol, perFreqID);
 
 // Get reference Epoch to time of first row or given String
-   
+
    Unit DAY(String("d"));
    MEpoch::Ref epochRef(in.getTimeReference());
@@ -1518 +1521 @@
      refEpoch = in.getEpoch(0);
    }
-   pushLog("Aligning at reference Epoch "+formatEpoch(refEpoch) 
-                  +" in frame "+MFrequency::showType(freqSystem));
-   
+   ostringstream oss;
+   oss << "Aligned at reference Epoch " << formatEpoch(refEpoch) << " in frame " << MFrequency::showType(freqSystem);
+   pushLog(String(oss));
+
    // Get Reference Position
-   
+
    MPosition refPos = in.getAntennaPosition();
-   
+
    // Create FrequencyAligner Block. One FA for each possible
    // source/freqID (perFreqID=True) or source/IF (perFreqID=False)
    // combination
-   
+
    PtrBlock<FrequencyAligner<Float>* > a(dDesc.length());
-   generateFrequencyAligners (a, dDesc, in, nChan, freqSystem, refPos, 
+   generateFrequencyAligners(a, dDesc, in, nChan, freqSystem, refPos,
                               refEpoch, perFreqID);
-   
+
    // Generate and fill output Frequency Table.  WHen perFreqID=True,
    // there is one output FreqID for each entry in the SDDataDesc
    // table.  However, in perFreqID=False mode, there may be some
    // degeneracy, so we need a little translation map
-   
+
    SDFrequencyTable freqTabOut = in.getSDFreqTable();
    freqTabOut.setLength(0);
@@ -1548 +1552 @@
 
      // Get Aligned SC in Hz
-     
+
      SpectralCoordinate sC = a[i]->alignedSpectralCoordinate(linear);
      sC.setWorldAxisUnits(units);
-     
+
      // Add FreqID
-     
+
      uInt idx = freqTabOut.addFrequency(sC.referencePixel()[0],
-                                        sC.referenceValue()[0],
-                                        sC.increment()[0]);
+                                        sC.referenceValue()[0],
+                                        sC.increment()[0]);
      // output FreqID = ddFQTrans(ddIdx)
      ddFQTrans(i) = idx;
    }
-   
+
    // Interpolation method
-   
+
    InterpolateArray1D<Double,Float>::InterpolationMethod interp;
    convertInterpString(interp, methodStr);
-   
+
    // New output Table
-   
-   pushLog("Create output table");
+
+   //pushLog("Create output table");
    SDMemTable* pTabOut = new SDMemTable(in,True);
    pTabOut->putSDFreqTable(freqTabOut);
-   
+
    // Loop over rows in Table
-   
+
    Bool extrapolate=False;
    const IPosition polChanAxes(2, asap::PolAxis, asap::ChanAxis);
@@ -1586 +1590 @@
    //
    for (uInt iRow=0; iRow<nRows; ++iRow) {
-     if (iRow%10==0) {
-       pushLog("Processing row "+iRow);
-     }
-     
+//      if (iRow%10==0) {
+//        ostringstream oss;
+//        oss << "Processing row " << iRow;
+//        pushLog(String(oss));
+//      }
+
      // Get EPoch
-     
+
      Quantum<Double> tQ2(times[iRow],DAY);
      MVEpoch mv2(tQ2);
      MEpoch epoch(mv2, epochRef);
-     
+
      // Get copy of data
-     
+
      const MaskedArray<Float>& mArrIn(in.rowAsMaskedArray(iRow));
      Array<Float> values = mArrIn.getArray();
-     Array<Bool> mask = mArrIn.getMask(); 
-     
+     Array<Bool> mask = mArrIn.getMask();
+
      // For each row, the Frequency abcissa will be the same
      // regardless of polarization.  For all other axes (IF and BEAM)
@@ -1607 +1613 @@
      // pol-chan planes to mimimize the work.  Probably won't work for
      // multiple beams at this point.
-     
+
      ArrayIterator<Float> itValuesPlane(values, polChanAxes);
      ArrayIterator<Bool> itMaskPlane(mask, polChanAxes);
@@ -1613 +1619 @@
 
        // Find the IF index and then the FA PtrBlock index
-       
+
        const IPosition& pos = itValuesPlane.pos();
        ifIdx = pos(asap::IFAxis);
        faIdx = ddIdx(iRow,ifIdx);
-       
+
        // Generate abcissa for perIF.  Could cache this in a Matrix on
        // a per scan basis.  Pretty expensive doing it for every row.
-       
+
        if (!perFreqID) {
-         xIn.resize(nChan);
-         uInt fqID = dDesc.secID(ddIdx(iRow,ifIdx));
-         SpectralCoordinate sC = in.getSpectralCoordinate(fqID);
+         xIn.resize(nChan);
+         uInt fqID = dDesc.secID(ddIdx(iRow,ifIdx));
+         SpectralCoordinate sC = in.getSpectralCoordinate(fqID);
            Double w;
            for (uInt i=0; i<nChan; i++) {
-             sC.toWorld(w,Double(i));
+             sC.toWorld(w,Double(i));
               xIn[i] = w;
            }
        }
-       
-       VectorIterator<Float> itValuesVec(itValuesPlane.array(), 1); 
+
+       VectorIterator<Float> itValuesVec(itValuesPlane.array(), 1);
        VectorIterator<Bool> itMaskVec(itMaskPlane.array(), 1);
 
        // Iterate through the plane by vector and align
-       
+
         first = True;
         useCachedAbcissa=False;
-        while (!itValuesVec.pastEnd()) {     
-          if (perFreqID) {
-            ok = a[faIdx]->align (yOut, maskOut, itValuesVec.vector(),
-                                  itMaskVec.vector(), epoch, useCachedAbcissa,
-                                  interp, extrapolate); 
-          } else {
-            ok = a[faIdx]->align (yOut, maskOut, xIn, itValuesVec.vector(),
-                                  itMaskVec.vector(), epoch, useCachedAbcissa,
-                                  interp, extrapolate); 
-          }
-          //
-          itValuesVec.vector() = yOut;
-          itMaskVec.vector() = maskOut;
-          //
-          itValuesVec.next();
-          itMaskVec.next();
-          //
-          if (first) {
-            useCachedAbcissa = True;
-            first = False;
-          }
+        while (!itValuesVec.pastEnd()) {
+          if (perFreqID) {
+            ok = a[faIdx]->align (yOut, maskOut, itValuesVec.vector(),
+                                  itMaskVec.vector(), epoch, useCachedAbcissa,
+                                  interp, extrapolate);
+          } else {
+            ok = a[faIdx]->align (yOut, maskOut, xIn, itValuesVec.vector(),
+                                  itMaskVec.vector(), epoch, useCachedAbcissa,
+                                  interp, extrapolate);
+          }
+          //
+          itValuesVec.vector() = yOut;
+          itMaskVec.vector() = maskOut;
+          //
+          itValuesVec.next();
+          itMaskVec.next();
+          //
+          if (first) {
+            useCachedAbcissa = True;
+            first = False;
+          }
         }
-        //
-        itValuesPlane.next();
-        itMaskPlane.next();
+        //
+        itValuesPlane.next();
+        itMaskPlane.next();
      }
-     
+
      // Create SDContainer and put back
-     
+
      SDContainer sc = in.getSDContainer(iRow);
      putDataInSDC(sc, values, mask);
-     
+
      // Set output FreqIDs
-     
+
      for (uInt i=0; i<nIF; i++) {
        uInt idx = ddIdx(iRow,i);               // Index into SDDataDesc table
@@ -1681 +1687 @@
      pTabOut->putSDContainer(sc);
    }
-   
+
    // Now we must set the base and extra frames to the input frame
    std::vector<string> info = pTabOut->getCoordInfo();
@@ -1688 +1694 @@
    pTabOut->setCoordInfo(info);
 
-   // Clean up PointerBlock   
+   // Clean up PointerBlock
    for (uInt i=0; i<a.nelements(); i++) delete a[i];
 
@@ -1704 +1710 @@
   const Table& tabIn = in.table();
 
-  // Shape of input and output data 
+  // Shape of input and output data
   const IPosition& shapeIn = in.rowAsMaskedArray(0).shape();
 
@@ -1750 +1756 @@
 
 void SDMath::fillSDC(SDContainer& sc,
-                     const Array<Bool>& mask,
-                     const Array<Float>& data,
-                     const Array<Float>& tSys,
-                     Int scanID, Double timeStamp,
-                     Double interval, const String& sourceName,
-                     const Vector<uInt>& freqID)
+                     const Array<Bool>& mask,
+                     const Array<Float>& data,
+                     const Array<Float>& tSys,
+                     Int scanID, Double timeStamp,
+                     Double interval, const String& sourceName,
+                     const Vector<uInt>& freqID)
 {
 // Data and mask
@@ -1776 +1782 @@
 
 void SDMath::accumulate(Double& timeSum, Double& intSum, Int& nAccum,
-                        MaskedArray<Float>& sum, Array<Float>& sumSq, 
-                        Array<Float>& nPts, Array<Float>& tSysSum, 
+                        MaskedArray<Float>& sum, Array<Float>& sumSq,
+                        Array<Float>& nPts, Array<Float>& tSysSum,
                         Array<Float>& tSysSqSum,
-                        const Array<Float>& tSys, const Array<Float>& nInc, 
-                        const Vector<Bool>& mask, Double time, Double interval,
-                        const std::vector<CountedPtr<SDMemTable> >& in,
-                        uInt iTab, uInt iRow, uInt axis, 
-                        uInt nAxesSub, Bool useMask,
-                        WeightType wtType)
+                        const Array<Float>& tSys, const Array<Float>& nInc,
+                        const Vector<Bool>& mask, Double time, Double interval,
+                        const std::vector<CountedPtr<SDMemTable> >& in,
+                        uInt iTab, uInt iRow, uInt axis,
+                        uInt nAxesSub, Bool useMask,
+                        WeightType wtType)
 {
 
@@ -1810 +1816 @@
       ReadOnlyVectorIterator<Bool> itMask(maskIn, axis);
       Float fac = 1.0;
-      IPosition pos(nAxesSub,0);  
+      IPosition pos(nAxesSub,0);
 //
       while (!itData.pastEnd()) {
@@ -1846 +1852 @@
 // We will take a short cut and just use the value from the first channel
 // for now.
-// 
+//
       VectorIterator<Float> itData(valuesIn, axis);
       ReadOnlyVectorIterator<Float> itTSys(tSys, axis);
-      IPosition pos(nAxesSub,0);  
+      IPosition pos(nAxesSub,0);
 //
       Float fac = 1.0;
@@ -1923 +1929 @@
       }
    } else if (wtType==TSYS || wtType==TINTSYS) {
-   
+
 // Normalize each spectrum by sum(1/Tsys**2) (TSYS) or
 // sum(Tint/Tsys**2) (TINTSYS) where the pseudo
@@ -2000 +2006 @@
 
 
-void SDMath::convertInterpString(casa::InterpolateArray1D<Double,Float>::InterpolationMethod& type,  
+void SDMath::convertInterpString(casa::InterpolateArray1D<Double,Float>::InterpolationMethod& type,
                                  const casa::String& interp)
 {
@@ -2019 +2025 @@
 
 void SDMath::putDataInSDC(SDContainer& sc, const Array<Float>& data,
-                          const Array<Bool>& mask)
+                          const Array<Bool>& mask)
 {
     sc.putSpectrum(data);
@@ -2051 +2057 @@
 }
 
-void SDMath::scaleFromTable(SDMemTable* pTabOut, const SDMemTable& in, 
-                            const Table& tTable, const String& col0, 
+void SDMath::scaleFromTable(SDMemTable* pTabOut, const SDMemTable& in,
+                            const Table& tTable, const String& col0,
                             const String& col1,
                             const String& methodStr, Bool doAll,
@@ -2074 +2080 @@
    Vector<Float> yOut;
    Vector<Bool> maskOut;
-   InterpolateArray1D<Float,Float>::interpolate(yOut, maskOut, xOut, 
+   InterpolateArray1D<Float,Float>::interpolate(yOut, maskOut, xOut,
                                                 xIn, yIn, maskIn, intMethod,
                                                 True, True);
-// Apply 
+// Apply
 
    scaleByVector(pTabOut, in, doAll, Float(1.0)/yOut, doTsys);
@@ -2083 +2089 @@
 
 
-void SDMath::scaleByVector(SDMemTable* pTabOut, const SDMemTable& in, 
+void SDMath::scaleByVector(SDMemTable* pTabOut, const SDMemTable& in,
                            Bool doAll, const Vector<Float>& factor,
                            Bool doTSys)
@@ -2100 +2106 @@
 
 // Loop over rows and apply correction factor
-  
+
   const uInt axis = asap::ChanAxis;
   for (uInt i=0; i < in.nRow(); ++i) {
@@ -2193 +2199 @@
 
 
-void SDMath::generateFrequencyAligners(PtrBlock<FrequencyAligner<Float>* >& a, 
-                                       const SDDataDesc& dDesc,
-                                       const SDMemTable& in, uInt nChan,
-                                       MFrequency::Types system,
-                                       const MPosition& refPos,
-                                       const MEpoch& refEpoch,
-                                       Bool perFreqID)
+void SDMath::generateFrequencyAligners(PtrBlock<FrequencyAligner<Float>* >& a,
+                                       const SDDataDesc& dDesc,
+                                       const SDMemTable& in, uInt nChan,
+                                       MFrequency::Types system,
+                                       const MPosition& refPos,
+                                       const MEpoch& refEpoch,
+                                       Bool perFreqID)
 {
    for (uInt i=0; i<dDesc.length(); i++) {
@@ -2208 +2214 @@
       if (perFreqID) {
 
-// One aligner per source/FreqID pair.  
+// One aligner per source/FreqID pair.
 
          SpectralCoordinate sC = in.getSpectralCoordinate(ID);
@@ -2230 +2236 @@
    return range;
 }
-   
+
 
 Bool SDMath::rowInRange(uInt i, const Vector<uInt>& range) const

branches/Release12/src/SDMemTable.cc

@@ -1703 +1703 @@
      }
      ostringstream oss;
-     oss << "Replacing rest frequencies, one per IF, with given list : " << restFreqs;
+     oss << "Replaced rest frequencies, one per IF, with given list : " << restFreqs;
      sdft.deleteRestFrequencies();
      for (uInt i=0; i<nRestFreqs; i++) {
@@ -1715 +1715 @@
       idx = sdft.addRestFrequency(rf.getValue("Hz"));
       if (whichIF>=0) {
-         oss << "Selecting given rest frequency (" << restFreqs[0] << ") for IF " << whichIF << endl;
+         oss << "Selected given rest frequency (" << restFreqs[0] << ") for IF " << whichIF << endl;
       } else {
-         oss << "Selecting given rest frequency (" << restFreqs[0] << ") for all IFs" << endl;
+         oss << "Selected given rest frequency (" << restFreqs[0] << ") for all IFs" << endl;
       }
    }
@@ -1839 +1839 @@
 {
   MPosition mp = getAntennaPosition();
-
+  ostringstream oss;
+  oss << "Computed azimuth/elevation using " << endl
+      << mp << endl;
   for (uInt i=0; i<nRow();++i) {
     MEpoch me = getEpoch(i);
     MDirection md = getDirection(i,False);
+    oss  << " Time: " << getTime(i,False) << " Direction: " << formatDirection(md)
+         << endl << "     => ";
     MeasFrame frame(mp, me);
     Vector<Double> azel =
@@ -1850 +1854 @@
     azCol_.put(i,azel[0]);
     elCol_.put(i,azel[1]);
-  }
-}
+    oss << "azel: " << azel[0]/C::pi*180.0 << " "
+        << azel[1]/C::pi*180.0 << " (deg)" << endl;
+  }
+  pushLog(String(oss));
+}