Changeset 470

Timestamp:

02/18/05 11:11:53 (20 years ago)

Author:

kil064

Message:

handle STokes output now in SDFITS/MS formats. Getting messy.

File:

• trunk/src/SDWriter.cc (modified) (10 diffs)

trunk/src/SDWriter.cc

@@ -1 @@
- //#---------------------------------------------------------------------------
+//#---------------------------------------------------------------------------
 //# SDWriter.cc: ASAP class to write out single dish spectra.
 //#---------------------------------------------------------------------------
@@ -33 +33 @@
 
 #include <casa/aips.h>
-#include <casa/Arrays.h>
+#include <casa/Arrays/Array.h>
+#include <casa/Arrays/Vector.h>
 #include <casa/BasicSL/Complex.h>
 #include <casa/Utilities/CountedPtr.h>
+#include <casa/Utilities/Assert.h>
 
 #include <atnf/PKSIO/PKSMS2writer.h>
@@ -41 +43 @@
 
 #include <tables/Tables/ArrayColumn.h>
+
+
 
 
@@ -48 +52 @@
 #include "SDFITSImageWriter.h"
 #include "SDAsciiWriter.h"
+#include "SDPol.h"
 
 using namespace casa;
@@ -147 +152 @@
   const Int nPol  = hdr.npol;
   const Int nChan = hdr.nchan;
-  if (toStokes) {
-     cerr << "Stokes conversion not yet available with SDFITS or MS" << endl;
-  }
 
 // Get Freq table
@@ -182 +184 @@
   Int cycleNo;
   Double mjd0 = 0.0;
-
+//
+  Array<Float> spectra, tSys, stokes;
+  Array<uChar> flags;
+  Bool doLinear = True;
+//
   Int count = 0;
   for (Int iRow = 0; iRow < in->nRow(); iRow++) {
@@ -198 +204 @@
     freqIDCol.get(iRow, freqIDs);
 
+    // Get Stokes array (all beams/IFs/Stokes).  Its not in SDContainer
+    if (toStokes) {
+       stokes = in->getStokesSpectrum(iRow,-1,-1);
+    }
+
+    IPosition start(asap::nAxes,0);
+    IPosition end(stokes.shape()-1);
+
     // Write it out beam by beam.
     for (Int iBeam = 0; iBeam < hdr.nbeam; iBeam++) {
+       start(asap::BeamAxis) = iBeam;
+       end(asap::BeamAxis) = iBeam;
 
       // Write it out IF by IF.
       for (Int iIF = 0; iIF < hdr.nif; iIF++) {
+        start(asap::IFAxis) = iIF;
+        end(asap::IFAxis) = iIF;
         uInt freqID = freqIDs(iIF);
 
@@ -242 +260 @@
         Complex         xCalFctr;
         Vector<Complex> xPol;
+
+// Get data.  
+
+        if (toStokes) {
+
+// Big pain in the ass because
+//  1) Stokes vector may have less polarizations than input
+//  2) Stokes column virtual and not in SDContainer
+//  3) Tsys and FLags already selected for IF and Beam 
+//  3) Have to flip order
+
+           Array<Float> t = sd.getTsys(iBeam,iIF);
+           tSys = SDPolUtil::computeStokesTSysForWriter (t, doLinear);
+           Array<uChar> t2 = sd.getFlags(iBeam,iIF);
+           flags = SDPolUtil::computeStokesFlagsForWriter (t2, doLinear);
+           spectra = SDPolUtil::extractStokesForWriter (stokes,start,end);
+        } else {
+           tSys = sd.getTsys(iBeam,iIF);
+           flags = sd.getFlags(iBeam,iIF);
+           spectra = sd.getSpectrum(iBeam,iIF);        
+        }
 //
         if (status = cWriter->write(sd.scanid, cycleNo, sd.timestamp,
@@ -255 +294 @@
                                     sd.getDirection(iBeam),
                                     scanRate,
-                                    sd.getTsys(iBeam, iIF), sigma, calFctr,
+                                    tSys, sigma, calFctr,
                                     baseLin, baseSub,
-                                    sd.getSpectrum(iBeam, iIF),
-                                    sd.getFlags(iBeam, iIF),
+                                    spectra, flags,
                                     xCalFctr, xPol)) {
           cerr << "Error writing output file." << endl;
@@ -274 +312 @@
   return 0;
 }
+
+