Ignore:
Timestamp:
11/19/08 20:41:16 (16 years ago)
Author:
Malte Marquarding
Message:

update from livedata CVS

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/external/atnf/PKSIO/PKSMS2reader.cc

    r1427 r1452  
    2626//#                        Charlottesville, VA 22903-2475 USA
    2727//#
    28 //# $Id: PKSMS2reader.cc,v 19.13 2008-06-26 02:08:02 cal103 Exp $
     28//# $Id: PKSMS2reader.cc,v 19.21 2008-11-17 06:55:18 cal103 Exp $
    2929//#---------------------------------------------------------------------------
    3030//# Original: 2000/08/03, Mark Calabretta, ATNF
    3131//#---------------------------------------------------------------------------
    3232
    33 
    34 // AIPS++ includes.
     33#include <atnf/pks/pks_maths.h>
     34#include <atnf/PKSIO/PKSmsg.h>
     35#include <atnf/PKSIO/PKSrecord.h>
     36#include <atnf/PKSIO/PKSMS2reader.h>
     37
    3538#include <casa/stdio.h>
    3639#include <casa/Arrays/ArrayMath.h>
     
    3942#include <tables/Tables.h>
    4043
    41 // Parkes includes.
    42 #include <atnf/pks/pks_maths.h>
    43 #include <atnf/PKSIO/PKSMS2reader.h>
    44 
    45 
    4644//------------------------------------------------- PKSMS2reader::PKSMS2reader
    4745
     
    5149{
    5250  cMSopen = False;
     51
     52  // By default, messages are written to stderr.
     53  initMsg();
    5354}
    5455
     
    353354        const Bool getSpectra,
    354355        const Bool getXPol,
    355         const Bool getFeedPos)
     356        const Int  coordSys)
    356357{
    357358  if (!cMSopen) {
     
    432433  cGetXPol = cGetXPol && getXPol;
    433434
    434   // Get feed positions?  (Not available.)
    435   cGetFeedPos = False;
     435  // Coordinate system?  (Only equatorial available.)
     436  cCoordSys = 0;
    436437
    437438  return maxNChan;
     
    486487// Read the next data record.
    487488
    488 Int PKSMS2reader::read(MBrecord &MBrec)
     489Int PKSMS2reader::read(PKSrecord &pksrec)
    489490{
    490491  if (!cMSopen) {
     
    514515
    515516  // Renumerate scan no. Here still is 1-based
    516   MBrec.scanNo = cScanNoCol(cIdx) - cScanNoCol(0) + 1;
    517 
    518   if (MBrec.scanNo != cScanNo) {
     517  pksrec.scanNo = cScanNoCol(cIdx) - cScanNoCol(0) + 1;
     518
     519  if (pksrec.scanNo != cScanNo) {
    519520    // Start of new scan.
    520     cScanNo  = MBrec.scanNo;
     521    cScanNo  = pksrec.scanNo;
    521522    cCycleNo = 1;
    522523    cTime    = cTimeCol(cIdx);
     
    524525
    525526  Double time = cTimeCol(cIdx);
    526   MBrec.mjd      = time/86400.0;
    527   MBrec.interval = cIntervalCol(cIdx);
     527  pksrec.mjd      = time/86400.0;
     528  pksrec.interval = cIntervalCol(cIdx);
    528529
    529530  // Reconstruct the integration cycle number; due to small latencies the
    530531  // integration time is usually slightly less than the time between cycles,
    531532  // resetting cTime will prevent the difference from accumulating.
    532   cCycleNo += nint((time - cTime)/MBrec.interval);
    533   MBrec.cycleNo = cCycleNo;
    534   cTime   = time;
     533  cCycleNo += nint((time - cTime)/pksrec.interval);
     534  pksrec.cycleNo = cCycleNo;
     535  cTime = time;
    535536
    536537  Int fieldId = cFieldIdCol(cIdx);
    537   MBrec.fieldName = cFieldNameCol(fieldId);
     538  pksrec.fieldName = cFieldNameCol(fieldId);
    538539
    539540  Int srcId = cSrcIdCol(fieldId);
    540   MBrec.srcName = cSrcNameCol(srcId);
    541   MBrec.srcDir  = cSrcDirCol(srcId);
    542   MBrec.srcPM   = cSrcPMCol(srcId);
     541  pksrec.srcName = cSrcNameCol(srcId);
     542  pksrec.srcDir  = cSrcDirCol(srcId);
     543  pksrec.srcPM   = cSrcPMCol(srcId);
    543544
    544545  // Systemic velocity.
    545546  if (!cHaveSrcVel) {
    546     MBrec.srcVel = 0.0f;
     547    pksrec.srcVel = 0.0f;
    547548  } else {
    548     MBrec.srcVel = cSrcVelCol(srcId)(IPosition(1,0));
     549    pksrec.srcVel = cSrcVelCol(srcId)(IPosition(1,0));
    549550  }
    550551
    551552  // Observation type.
    552553  Int stateId = cStateIdCol(cIdx);
    553   MBrec.obsType = cObsModeCol(stateId);
    554 
    555   MBrec.IFno = iIF + 1;
     554  pksrec.obsType = cObsModeCol(stateId);
     555
     556  pksrec.IFno = iIF + 1;
    556557  Int nChan = abs(cEndChan(iIF) - cStartChan(iIF)) + 1;
    557558
     
    560561  if (nChan == 1) {
    561562    cout << "The input is continuum data. "<< endl;
    562     MBrec.freqInc  = chanFreq(0);
    563     MBrec.refFreq  = chanFreq(0);
    564     MBrec.restFreq = 0.0f;
     563    pksrec.freqInc  = chanFreq(0);
     564    pksrec.refFreq  = chanFreq(0);
     565    pksrec.restFreq = 0.0f;
    565566  } else {
    566567    if (cStartChan(iIF) <= cEndChan(iIF)) {
    567       MBrec.freqInc = chanFreq(1) - chanFreq(0);
     568      pksrec.freqInc = chanFreq(1) - chanFreq(0);
    568569    } else {
    569       MBrec.freqInc = chanFreq(0) - chanFreq(1);
    570     }
    571 
    572     MBrec.refFreq  = chanFreq(cRefChan(iIF)-1);
    573     MBrec.restFreq = cSrcRestFrqCol(srcId)(IPosition(1,0));
    574   }
    575   MBrec.bandwidth = abs(MBrec.freqInc * nChan);
    576 
    577   MBrec.tcal.resize(cNPol(iIF));
    578   MBrec.tcal      = 0.0f;
    579   MBrec.tcalTime  = "";
    580   MBrec.azimuth   = 0.0f;
    581   MBrec.elevation = 0.0f;
    582   MBrec.parAngle  = 0.0f;
    583   MBrec.focusAxi  = 0.0f;
    584   MBrec.focusTan  = 0.0f;
    585   MBrec.focusRot  = 0.0f;
     570      pksrec.freqInc = chanFreq(0) - chanFreq(1);
     571    }
     572
     573    pksrec.refFreq  = chanFreq(cRefChan(iIF)-1);
     574    pksrec.restFreq = cSrcRestFrqCol(srcId)(IPosition(1,0));
     575  }
     576  pksrec.bandwidth = abs(pksrec.freqInc * nChan);
     577
     578  pksrec.tcal.resize(cNPol(iIF));
     579  pksrec.tcal      = 0.0f;
     580  pksrec.tcalTime  = "";
     581  pksrec.azimuth   = 0.0f;
     582  pksrec.elevation = 0.0f;
     583  pksrec.parAngle  = 0.0f;
     584
     585  pksrec.focusAxi  = 0.0f;
     586  pksrec.focusTan  = 0.0f;
     587  pksrec.focusRot  = 0.0f;
    586588
    587589  // Find the appropriate entry in the WEATHER subtable.
     
    596598  if (weatherIdx < 0) {
    597599    // No appropriate WEATHER entry.
    598     MBrec.temperature = 0.0f;
    599     MBrec.pressure    = 0.0f;
    600     MBrec.humidity    = 0.0f;
     600    pksrec.temperature = 0.0f;
     601    pksrec.pressure    = 0.0f;
     602    pksrec.humidity    = 0.0f;
    601603  } else {
    602     MBrec.temperature = cTemperatureCol(weatherIdx);
    603     MBrec.pressure    = cPressureCol(weatherIdx);
    604     MBrec.humidity    = cHumidityCol(weatherIdx);
    605   }
    606 
    607   MBrec.windSpeed = 0.0f;
    608   MBrec.windAz    = 0.0f;
    609 
    610   MBrec.refBeam = 0;
    611   MBrec.beamNo  = ibeam + 1;
     604    pksrec.temperature = cTemperatureCol(weatherIdx);
     605    pksrec.pressure    = cPressureCol(weatherIdx);
     606    pksrec.humidity    = cHumidityCol(weatherIdx);
     607  }
     608
     609  pksrec.windSpeed = 0.0f;
     610  pksrec.windAz    = 0.0f;
     611
     612  pksrec.refBeam = 0;
     613  pksrec.beamNo  = ibeam + 1;
    612614
    613615  Matrix<Double> pointingDir = cPointingCol(fieldId);
    614   MBrec.direction = pointingDir.column(0);
     616  pksrec.direction = pointingDir.column(0);
     617  pksrec.pCode = 0;
     618  pksrec.rateAge = 0.0f;
    615619  uInt ncols = pointingDir.ncolumn();
    616620  if (ncols == 1) {
    617     MBrec.scanRate = 0.0f;
     621    pksrec.scanRate = 0.0f;
    618622  } else {
    619     MBrec.scanRate  = pointingDir.column(1);
    620   }
    621   MBrec.rateAge = 0;
    622   MBrec.rateson = 0;
     623    pksrec.scanRate(0) = pointingDir.column(1)(0);
     624    pksrec.scanRate(1) = pointingDir.column(1)(1);
     625  }
     626  pksrec.paRate = 0.0f;
    623627
    624628  // Get Tsys assuming that entries in the SYSCAL table match the main table.
     
    630634  }
    631635  if (cHaveTsys) {
    632     cTsysCol.get(cIdx, MBrec.tsys, True);
     636    cTsysCol.get(cIdx, pksrec.tsys, True);
    633637  } else {
    634638    Int numReceptor;
    635639    cNumReceptorCol.get(0, numReceptor);
    636     MBrec.tsys.resize(numReceptor);
    637     MBrec.tsys = 1.0f;
    638   }
    639   cSigmaCol.get(cIdx, MBrec.sigma, True);
     640    pksrec.tsys.resize(numReceptor);
     641    pksrec.tsys = 1.0f;
     642  }
     643  cSigmaCol.get(cIdx, pksrec.sigma, True);
    640644
    641645  // Calibration factors (if available).
    642   MBrec.calFctr.resize(cNPol(iIF));
     646  pksrec.calFctr.resize(cNPol(iIF));
    643647  if (cHaveCalFctr) {
    644     cCalFctrCol.get(cIdx, MBrec.calFctr);
     648    cCalFctrCol.get(cIdx, pksrec.calFctr);
    645649  } else {
    646     MBrec.calFctr = 0.0f;
     650    pksrec.calFctr = 0.0f;
    647651  }
    648652
    649653  // Baseline parameters (if available).
    650654  if (cHaveBaseLin) {
    651     MBrec.baseLin.resize(2,cNPol(iIF));
    652     cBaseLinCol.get(cIdx, MBrec.baseLin);
    653 
    654     MBrec.baseSub.resize(9,cNPol(iIF));
    655     cBaseSubCol.get(cIdx, MBrec.baseSub);
     655    pksrec.baseLin.resize(2,cNPol(iIF));
     656    cBaseLinCol.get(cIdx, pksrec.baseLin);
     657
     658    pksrec.baseSub.resize(9,cNPol(iIF));
     659    cBaseSubCol.get(cIdx, pksrec.baseSub);
    656660
    657661  } else {
    658     MBrec.baseLin.resize(0,0);
    659     MBrec.baseSub.resize(0,0);
     662    pksrec.baseLin.resize(0,0);
     663    pksrec.baseSub.resize(0,0);
    660664  }
    661665
     
    670674    // Transpose spectra.
    671675    Int nPol = tmpData.nrow();
    672     MBrec.spectra.resize(nChan, nPol);
    673     MBrec.flagged.resize(nChan, nPol);
     676    pksrec.spectra.resize(nChan, nPol);
     677    pksrec.flagged.resize(nChan, nPol);
    674678    if (cEndChan(iIF) >= cStartChan(iIF)) {
    675679      // Simple transposition.
    676680      for (Int ipol = 0; ipol < nPol; ipol++) {
    677681        for (Int ichan = 0; ichan < nChan; ichan++) {
    678           MBrec.spectra(ichan,ipol) = tmpData(ipol,ichan);
    679           MBrec.flagged(ichan,ipol) = tmpFlag(ipol,ichan);
     682          pksrec.spectra(ichan,ipol) = tmpData(ipol,ichan);
     683          pksrec.flagged(ichan,ipol) = tmpFlag(ipol,ichan);
    680684        }
    681685      }
     
    686690      for (Int ipol = 0; ipol < nPol; ipol++) {
    687691        for (Int ichan = 0; ichan < nChan; ichan++, jchan--) {
    688           MBrec.spectra(ichan,ipol) = tmpData(ipol,jchan);
    689           MBrec.flagged(ichan,ipol) = tmpFlag(ipol,jchan);
     692          pksrec.spectra(ichan,ipol) = tmpData(ipol,jchan);
     693          pksrec.flagged(ichan,ipol) = tmpFlag(ipol,jchan);
    690694        }
    691695      }
     
    696700  if (cGetXPol) {
    697701    if (cHaveXCalFctr) {
    698       cXCalFctrCol.get(cIdx, MBrec.xCalFctr);
     702      cXCalFctrCol.get(cIdx, pksrec.xCalFctr);
    699703    } else {
    700       MBrec.xCalFctr = Complex(0.0f, 0.0f);
    701     }
    702 
    703     cDataCol.get(cIdx, MBrec.xPol, True);
     704      pksrec.xCalFctr = Complex(0.0f, 0.0f);
     705    }
     706
     707    cDataCol.get(cIdx, pksrec.xPol, True);
    704708
    705709    if (cEndChan(iIF) < cStartChan(iIF)) {
     
    707711      Int jchan = nChan - 1;
    708712      for (Int ichan = 0; ichan < nChan/2; ichan++, jchan--) {
    709         ctmp = MBrec.xPol(ichan);
    710         MBrec.xPol(ichan) = MBrec.xPol(jchan);
    711         MBrec.xPol(jchan) = ctmp;
    712       }
    713     }
    714   }
    715 
    716   cIdx++;
    717 
    718   return 0;
    719 }
    720 
    721 //--------------------------------------------------------- PKSMS2reader::read
    722 
    723 // Read the next data record, just the basics.
    724 
    725 Int PKSMS2reader::read(
    726         Int           &IFno,
    727         Vector<Float> &tsys,
    728         Vector<Float> &calFctr,
    729         Matrix<Float> &baseLin,
    730         Matrix<Float> &baseSub,
    731         Matrix<Float> &spectra,
    732         Matrix<uChar> &flagged)
    733 {
    734   if (!cMSopen) {
    735     return 1;
    736   }
    737 
    738   // Check for EOF.
    739   if (cIdx >= cNRow) {
    740     return -1;
    741   }
    742 
    743   // Find the next selected beam and IF.
    744   Int ibeam;
    745   Int iIF;
    746   while (True) {
    747     ibeam = cBeamNoCol(cIdx);
    748     iIF   = cDataDescIdCol(cIdx);
    749     if (cBeams(ibeam) && cIFs(iIF)) {
    750       break;
    751     }
    752 
    753     // Check for EOF.
    754     if (++cIdx >= cNRow) {
    755       return -1;
    756     }
    757   }
    758 
    759   IFno = iIF + 1;
    760 
    761   // Get Tsys assuming that entries in the SYSCAL table match the main table.
    762   cTsysCol.get(cIdx, tsys, True);
    763 
    764   // Calibration factors (if available).
    765   if (cHaveCalFctr) {
    766     cCalFctrCol.get(cIdx, calFctr, True);
    767   } else {
    768     calFctr.resize(cNPol(iIF));
    769     calFctr = 0.0f;
    770   }
    771 
    772   // Baseline parameters (if available).
    773   if (cHaveBaseLin) {
    774     baseLin.resize(2,cNPol(iIF));
    775     cBaseLinCol.get(cIdx, baseLin);
    776 
    777     baseSub.resize(9,cNPol(iIF));
    778     cBaseSubCol.get(cIdx, baseSub);
    779 
    780   } else {
    781     baseLin.resize(0,0);
    782     baseSub.resize(0,0);
    783   }
    784 
    785   if (cGetSpectra) {
    786     // Get spectral data.
    787     Matrix<Float> tmpData;
    788     Matrix<Bool>  tmpFlag;
    789     cFloatDataCol.getSlice(cIdx, cDataSel(iIF), tmpData, True);
    790     cFlagCol.getSlice(cIdx, cDataSel(iIF), tmpFlag, True);
    791 
    792     // Transpose spectra.
    793     Int nChan = tmpData.ncolumn();
    794     Int nPol  = tmpData.nrow();
    795     spectra.resize(nChan, nPol);
    796     flagged.resize(nChan, nPol);
    797     if (cEndChan(iIF) >= cStartChan(iIF)) {
    798       // Simple transposition.
    799       for (Int ipol = 0; ipol < nPol; ipol++) {
    800         for (Int ichan = 0; ichan < nChan; ichan++) {
    801           spectra(ichan,ipol) = tmpData(ipol,ichan);
    802           flagged(ichan,ipol) = tmpFlag(ipol,ichan);
    803         }
    804       }
    805 
    806     } else {
    807       // Transpose with inversion.
    808       Int jchan = nChan - 1;
    809       for (Int ipol = 0; ipol < nPol; ipol++) {
    810         for (Int ichan = 0; ichan < nChan; ichan++, jchan--) {
    811           spectra(ichan,ipol) = tmpData(ipol,jchan);
    812           flagged(ichan,ipol) = tmpFlag(ipol,jchan);
    813         }
     713        ctmp = pksrec.xPol(ichan);
     714        pksrec.xPol(ichan) = pksrec.xPol(jchan);
     715        pksrec.xPol(jchan) = ctmp;
    814716      }
    815717    }
Note: See TracChangeset for help on using the changeset viewer.