Ignore:
Timestamp:
07/18/09 06:35:47 (15 years ago)
Author:
TakTsutsumi
Message:

New Development: No, merge with asap2.3.1

JIRA Issue: Yes CAS-1450

Ready to Release: Yes/No?

Interface Changes: Yes/No?

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: Yes

Module(s): single dish

Description: Upgrade of alma branch based on ASAP2.2.0

(rev.1562) to ASAP2.3.1 (rev.1561)


File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/alma/src/STWriter.cpp

    r1446 r1603  
    3939#include <casa/Utilities/Assert.h>
    4040
     41#include <atnf/PKSIO/PKSrecord.h>
    4142#include <atnf/PKSIO/PKSMS2writer.h>
    4243#include <atnf/PKSIO/PKSSDwriter.h>
    … …  
    4748#include <tables/Tables/ArrayColumn.h>
    4849
    49 //#include "SDFITSImageWriter.h"
     50#include "STFITSImageWriter.h"
    5051#include "STAsciiWriter.h"
    5152#include "STHeader.h"
    … …  
    5859STWriter::STWriter(const std::string &format)
    5960{
     61  format_ = format;
     62  String t(format_);
     63  t.upcase();
     64  if (t == "MS2") {
     65    writer_ = new PKSMS2writer();
     66  } else if (t == "SDFITS") {
     67    writer_ = new PKSSDwriter();
     68  } else if (t == "ASCII" || t == "FITS" || t == "CLASS") {
     69    writer_ = 0;
     70  } else {
     71    throw (AipsError("Unrecognized export format"));
     72  }
     73}
     74
     75STWriter::~STWriter()
     76{
     77   if (writer_) {
     78     delete writer_;
     79   }
     80}
     81
     82Int STWriter::setFormat(const std::string &format)
     83{
     84  if (format != format_) {
     85    if (writer_) delete writer_;
     86  }
     87
    6088  format_ = format;
    6189  String t(format_);
    … …  
    6593  } else if (t== "SDFITS") {
    6694    writer_ = new PKSSDwriter();
    67   } else if (t== "ASCII") {
    68     writer_ = 0;
    69   } else {
    70     throw (AipsError("Unrecognized export format"));
    71   }
    72 }
    73 
    74 STWriter::~STWriter()
    75 {
    76    if (writer_) {
    77      delete writer_;
    78    }
    79 }
    80 
    81 Int STWriter::setFormat(const std::string &format)
    82 {
    83   if (format != format_) {
    84     if (writer_) delete writer_;
    85   }
    86 
    87   format_ = format;
    88   String t(format_);
    89   t.upcase();
    90   if (t== "MS2") {
    91     writer_ = new PKSMS2writer();
    92   } else if (t== "SDFITS") {
    93     writer_ = new PKSSDwriter();
    94   } else if (t== "ASCII") {
     95  } else if (t == "ASCII" || t == "FITS" || t == "CLASS") {
    9596    writer_ = 0;
    9697  } else {
    … …  
    110111    } else {
    111112      return 1;
     113    }
     114  } else if ( format_ == "FITS" || format_ == "CLASS") {
     115    STFITSImageWriter iw;
     116    if (format_ == "CLASS") {
     117      iw.setClass(True);
     118    }
     119    if (iw.write(*in, filename)) {
     120      return 0;
    112121    }
    113122  }
    … …  
    139148  Int status;
    140149  status = writer_->create(String(filename), hdr.observer, hdr.project,
    141                                hdr.antennaname, hdr.antennaposition,
    142                                hdr.obstype, hdr.equinox, hdr.freqref,
    143                                nChan, nPol, havexpol, False, fluxUnit);
     150                           hdr.antennaname, hdr.antennaposition,
     151                           hdr.obstype, hdr.fluxunit,
     152                           hdr.equinox, hdr.freqref,
     153                           nChan, nPol, havexpol, False);
    144154  if ( status ) {
    145155    throw(AipsError("Failed to create output file"));
    146156  }
    147157
    148   Double          srcVel;
    149 
    150   String          fieldName, srcName, tcalTime;
    151   Vector<Float>   calFctr, sigma, tcal, tsys;
    152   Vector<Double>  direction(2), scanRate(2), srcDir(2), srcPM(2);
    153   Matrix<Float>   spectra;
    154   Matrix<uChar>   flagtra;
    155   Complex         xCalFctr;
     158
    156159  Int count = 0;
    157   Int scanno = 1;
     160  PKSrecord pksrec;
     161  pksrec.scanNo = 1;
    158162  // use spearate iterators to ensure renumbering of all numbers
    159163  TableIterator scanit(table, "SCANNO");
    … …  
    161165    Table stable = scanit.table();
    162166    TableIterator beamit(stable, "BEAMNO");
    163     Int beamno = 1;
     167    pksrec.beamNo = 1;
    164168    while (!beamit.pastEnd() ) {
    165169      Table btable = beamit.table();
    166       // position only varies by beam
    167       // No, we want to pointing data which varies by cycle!
    168170      MDirection::ScalarColumn dirCol(btable, "DIRECTION");
    169       Vector<Double> direction = dirCol(0).getAngle("rad").getValue();
     171      pksrec.direction = dirCol(0).getAngle("rad").getValue();
    170172      TableIterator cycit(btable, "CYCLENO");
    171173      ROArrayColumn<Double> srateCol(btable, "SCANRATE");
    172       srateCol.get(0, scanRate);
     174      Vector<Double> sratedbl;
     175      srateCol.get(0, sratedbl);
     176      Vector<Float> srateflt(sratedbl.nelements());
     177      convertArray(srateflt, sratedbl);
     178      //pksrec.scanRate = srateflt;
     179      pksrec.scanRate = sratedbl;
    173180      ROArrayColumn<Double> spmCol(btable, "SRCPROPERMOTION");
    174       spmCol.get(0, srcPM);
     181      spmCol.get(0, pksrec.srcPM);
    175182      ROArrayColumn <Double> sdirCol(btable, "SRCDIRECTION");
    176       sdirCol.get(0, srcDir);
     183      sdirCol.get(0, pksrec.srcDir);
    177184      ROScalarColumn<Double> svelCol(btable, "SRCVELOCITY");
    178       svelCol.get(0, srcVel);
     185      svelCol.get(0, pksrec.srcVel);
    179186      ROScalarColumn<uInt> bCol(btable, "BEAMNO");
    180       beamno = bCol(0)+1;
    181       Int cycno = 1;
     187      pksrec.beamNo = bCol(0)+1;
     188      pksrec.cycleNo = 1;
    182189      while (!cycit.pastEnd() ) {
    183190        Table ctable = cycit.table();
     191        TableIterator ifit(ctable, "IFNO");
    184192        //MDirection::ScalarColumn dirCol(ctable, "DIRECTION");
    185         //Vector<Double> direction = dirCol(0).getAngle("rad").getValue();
    186         TableIterator ifit(ctable, "IFNO");
    187         Int ifno = 1;
     193        //pksrec.direction = dirCol(0).getAngle("rad").getValue();
     194        pksrec.IFno = 1;
    188195        while (!ifit.pastEnd() ) {
    189196          Table itable = ifit.table();
    … …  
    192199          const TableRecord& rec = row.get(0);
    193200          ROArrayColumn<Float> specCol(itable, "SPECTRA");
    194           ifno = rec.asuInt("IFNO")+1;
     201          pksrec.IFno = rec.asuInt("IFNO")+1;
    195202          uInt nchan = specCol(0).nelements();
    196           //Double cdelt,crval,crpix, restfreq;
    197           Double cdelt,crval,crpix;
    198           Vector<Double> restfreq;
    199           Float focusAxi, focusTan, focusRot,
    200                 temperature, pressure, humidity, windSpeed, windAz;
     203          Double crval,crpix;
     204          //Vector<Double> restfreq;
    201205          Float tmp0,tmp1,tmp2,tmp3,tmp4;
    202           Vector<Float> tcalval;
    203           //String stmp0,stmp1, tcalt;
    204206          String tcalt;
    205207          Vector<String> stmp0, stmp1;
    206           in->frequencies().getEntry(crpix,crval,cdelt, rec.asuInt("FREQ_ID"));
    207           in->focus().getEntry(focusAxi, focusTan, focusRot,
    208                                tmp0,tmp1,tmp2,tmp3,tmp4,
     208          in->frequencies().getEntry(crpix,crval, pksrec.freqInc,
     209                                     rec.asuInt("FREQ_ID"));
     210          in->focus().getEntry(pksrec.focusAxi, pksrec.focusTan,
     211                               pksrec.focusRot, tmp0,tmp1,tmp2,tmp3,tmp4,
    209212                               rec.asuInt("FOCUS_ID"));
    210           in->molecules().getEntry(restfreq,stmp0,stmp1,rec.asuInt("MOLECULE_ID"));
    211           in->tcal().getEntry(tcalt,tcalval,rec.asuInt("TCAL_ID"));
    212           in->weather().getEntry(temperature, pressure, humidity,
    213                                  windSpeed, windAz,
    214                                  rec.asuInt("WEATHER_ID"));
     213          in->molecules().getEntry(pksrec.restFreq,stmp0,stmp1,
     214                                   rec.asuInt("MOLECULE_ID"));
     215          in->tcal().getEntry(pksrec.tcalTime, pksrec.tcal,
     216                              rec.asuInt("TCAL_ID"));
     217          in->weather().getEntry(pksrec.temperature, pksrec.pressure,
     218                                 pksrec.humidity, pksrec.windSpeed,
     219                                 pksrec.windAz, rec.asuInt("WEATHER_ID"));
    215220          Double pixel = Double(nchan/2);
    216           Double refFreqNew = (pixel-crpix)*cdelt + crval;
     221          pksrec.refFreq = (pixel-crpix)*pksrec.freqInc + crval;
    217222          // ok, now we have nrows for the n polarizations in this table
    218           Matrix<Float> specs;
    219           Matrix<uChar> flags;
    220           Vector<Complex> xpol;
    221           polConversion(specs, flags, xpol, itable);
    222           Vector<Float> tsys = tsysFromTable(itable);
     223          polConversion(pksrec.spectra, pksrec.flagged, pksrec.xPol, itable);
     224          pksrec.tsys = tsysFromTable(itable);
    223225          // dummy data
    224           //uInt npol;
    225           //if ( hdr.antennaname == "GBT" ) {
    226           //  npol = nPolUsed;
    227           //}
    228           //else {   
    229           //  npol = specs.ncolumn();
    230           //}
    231           uInt npol = specs.ncolumn();
    232 
    233           Matrix<Float>   baseLin(npol,2, 0.0f);
    234           Matrix<Float>   baseSub(npol,9, 0.0f);
    235           Complex         xCalFctr;
    236           Vector<Double>  scanRate(2, 0.0);
    237           Vector<Float>   sigma(npol, 0.0f);
    238           Vector<Float>   calFctr(npol, 0.0f);
    239           status = writer_->write(scanno, cycno, rec.asDouble("TIME"),
    240                                       rec.asDouble("INTERVAL"),
    241                                       rec.asString("FIELDNAME"),
    242                                       rec.asString("SRCNAME"),
    243                                       srcDir, srcPM, srcVel,hdr.obstype,
    244                                       ifno,
    245                                       refFreqNew, nchan*abs(cdelt), cdelt,
    246                                       restfreq,
    247                                       tcalval,
    248                                       tcalt,
    249                                       rec.asFloat("AZIMUTH"),
    250                                       rec.asFloat("ELEVATION"),
    251                                       rec.asFloat("PARANGLE"),
    252                                       focusAxi, focusTan, focusRot,
    253                                       temperature,
    254                                       pressure, humidity, windSpeed, windAz,
    255                                       rec.asInt("REFBEAMNO")+1, beamno,
    256                                       direction,
    257                                       scanRate,
    258                                       tsys,
    259                                       sigma, calFctr,// not in scantable
    260                                       baseLin, baseSub,// not in scantable
    261                                       specs, flags,
    262                                       xCalFctr,//
    263                                       xpol);
     226          uInt npol = pksrec.spectra.ncolumn();
     227
     228          pksrec.mjd       = rec.asDouble("TIME");
     229          pksrec.interval  = rec.asDouble("INTERVAL");
     230          pksrec.fieldName = rec.asString("FIELDNAME");
     231          pksrec.srcName   = rec.asString("SRCNAME");
     232          pksrec.obsType   = hdr.obstype;
     233          pksrec.bandwidth = nchan * abs(pksrec.freqInc);
     234          pksrec.azimuth   = rec.asFloat("AZIMUTH");
     235          pksrec.elevation = rec.asFloat("ELEVATION");
     236          pksrec.parAngle  = rec.asFloat("PARANGLE");
     237          pksrec.refBeam   = rec.asInt("REFBEAMNO") + 1;
     238          pksrec.sigma.resize(npol);
     239          pksrec.sigma     = 0.0f;
     240          pksrec.calFctr.resize(npol);
     241          pksrec.calFctr   = 0.0f;
     242          pksrec.baseLin.resize(npol,2);
     243          pksrec.baseLin   = 0.0f;
     244          pksrec.baseSub.resize(npol,9);
     245          pksrec.baseSub   = 0.0f;
     246          pksrec.xCalFctr  = 0.0;
     247
     248          status = writer_->write(pksrec);
    264249          if ( status ) {
    265250            writer_->close();
    … …  
    267252          }
    268253          ++count;
    269           ++ifno;
     254          //++pksrec.IFno;
    270255          ++ifit;
    271256        }
    272         ++cycno;
     257        ++pksrec.cycleNo;
    273258        ++cycit;
    274259      }
    275       ++beamno;
     260      //++pksrec.beamNo;
    276261      ++beamit;
    277262    }
    278     ++scanno;
     263    ++pksrec.scanNo;
    279264    ++scanit;
    280265  }
    … …  
    286271  if ( format_ == "MS2" ) {
    287272    replacePtTab(table, filename);
    288   } 
     273  }
    289274  return 0;
    290275}
Note: See TracChangeset for help on using the changeset viewer.