Changeset 470 for trunk


Ignore:
Timestamp:
02/18/05 11:11:53 (20 years ago)
Author:
kil064
Message:

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

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/SDWriter.cc

    r450 r470  
    1  //#---------------------------------------------------------------------------
     1//#---------------------------------------------------------------------------
    22//# SDWriter.cc: ASAP class to write out single dish spectra.
    33//#---------------------------------------------------------------------------
     
    3333
    3434#include <casa/aips.h>
    35 #include <casa/Arrays.h>
     35#include <casa/Arrays/Array.h>
     36#include <casa/Arrays/Vector.h>
    3637#include <casa/BasicSL/Complex.h>
    3738#include <casa/Utilities/CountedPtr.h>
     39#include <casa/Utilities/Assert.h>
    3840
    3941#include <atnf/PKSIO/PKSMS2writer.h>
     
    4143
    4244#include <tables/Tables/ArrayColumn.h>
     45
     46
    4347
    4448
     
    4852#include "SDFITSImageWriter.h"
    4953#include "SDAsciiWriter.h"
     54#include "SDPol.h"
    5055
    5156using namespace casa;
     
    147152  const Int nPol  = hdr.npol;
    148153  const Int nChan = hdr.nchan;
    149   if (toStokes) {
    150      cerr << "Stokes conversion not yet available with SDFITS or MS" << endl;
    151   }
    152154
    153155// Get Freq table
     
    182184  Int cycleNo;
    183185  Double mjd0 = 0.0;
    184 
     186//
     187  Array<Float> spectra, tSys, stokes;
     188  Array<uChar> flags;
     189  Bool doLinear = True;
     190//
    185191  Int count = 0;
    186192  for (Int iRow = 0; iRow < in->nRow(); iRow++) {
     
    198204    freqIDCol.get(iRow, freqIDs);
    199205
     206    // Get Stokes array (all beams/IFs/Stokes).  Its not in SDContainer
     207    if (toStokes) {
     208       stokes = in->getStokesSpectrum(iRow,-1,-1);
     209    }
     210
     211    IPosition start(asap::nAxes,0);
     212    IPosition end(stokes.shape()-1);
     213
    200214    // Write it out beam by beam.
    201215    for (Int iBeam = 0; iBeam < hdr.nbeam; iBeam++) {
     216       start(asap::BeamAxis) = iBeam;
     217       end(asap::BeamAxis) = iBeam;
    202218
    203219      // Write it out IF by IF.
    204220      for (Int iIF = 0; iIF < hdr.nif; iIF++) {
     221        start(asap::IFAxis) = iIF;
     222        end(asap::IFAxis) = iIF;
    205223        uInt freqID = freqIDs(iIF);
    206224
     
    242260        Complex         xCalFctr;
    243261        Vector<Complex> xPol;
     262
     263// Get data. 
     264
     265        if (toStokes) {
     266
     267// Big pain in the ass because
     268//  1) Stokes vector may have less polarizations than input
     269//  2) Stokes column virtual and not in SDContainer
     270//  3) Tsys and FLags already selected for IF and Beam
     271//  3) Have to flip order
     272
     273           Array<Float> t = sd.getTsys(iBeam,iIF);
     274           tSys = SDPolUtil::computeStokesTSysForWriter (t, doLinear);
     275           Array<uChar> t2 = sd.getFlags(iBeam,iIF);
     276           flags = SDPolUtil::computeStokesFlagsForWriter (t2, doLinear);
     277           spectra = SDPolUtil::extractStokesForWriter (stokes,start,end);
     278        } else {
     279           tSys = sd.getTsys(iBeam,iIF);
     280           flags = sd.getFlags(iBeam,iIF);
     281           spectra = sd.getSpectrum(iBeam,iIF);       
     282        }
    244283//
    245284        if (status = cWriter->write(sd.scanid, cycleNo, sd.timestamp,
     
    255294                                    sd.getDirection(iBeam),
    256295                                    scanRate,
    257                                     sd.getTsys(iBeam, iIF), sigma, calFctr,
     296                                    tSys, sigma, calFctr,
    258297                                    baseLin, baseSub,
    259                                     sd.getSpectrum(iBeam, iIF),
    260                                     sd.getFlags(iBeam, iIF),
     298                                    spectra, flags,
    261299                                    xCalFctr, xPol)) {
    262300          cerr << "Error writing output file." << endl;
     
    274312  return 0;
    275313}
     314
     315
Note: See TracChangeset for help on using the changeset viewer.