Ignore:
Timestamp:
11/27/08 12:47:15 (16 years ago)
Author:
TakTsutsumi
Message:

New Development: No

JIRA Issue: No

Ready to Release: Yes

Interface Changes: Yes

What Interface Changed: many

Test Programs: sd.scantable(), sd.scantable.save()

Put in Release Notes: N/A

Description: copied from current casapy code tree


File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/alma/external/atnf/PKSIO/PKSMS2reader.cc

    r1393 r1453  
    9292  cPKSMS  = MeasurementSet(msName);
    9393  // taql access to the syscal table
    94   String tmpcalName = msName;
    95   tmpcalName.append("/SYSCAL");
    96   cSysCalTab = Table(tmpcalName);
     94  cHaveSysCal = False;
     95  if (cHaveSysCal=Table::isReadable(cPKSMS.sysCalTableName())) {
     96    cSysCalTab = Table(cPKSMS.sysCalTableName());
     97  }
    9798
    9899  cIdx    = 0;
    … …  
    116117  ROMSSysCalColumns       sysCalCols(cPKSMS.sysCal());
    117118  ROMSWeatherColumns      weatherCols(cPKSMS.weather());
     119  ROMSAntennaColumns      antennaCols(cPKSMS.antenna());
    118120
    119121  // Column accessors for required columns.
    … …  
    153155
    154156  // Optional columns.
     157  cHaveTsys = False;
     158  cHaveTcal = False;
    155159  if ((cHaveSrcVel = cPKSMS.source().tableDesc().isColumn("SYSVEL"))) {
    156160    cSrcVelCol.attach(cPKSMS.source(), "SYSVEL");
    157161  }
    158162
    159   if ((cHaveTsys = cPKSMS.sysCal().tableDesc().isColumn("TSYS"))) {
     163  if (cHaveSysCal && (cHaveTsys = cPKSMS.sysCal().tableDesc().isColumn("TSYS"))) {
    160164    cTsysCol.attach(cPKSMS.sysCal(), "TSYS");
    161165  }
    162166 
    163   if ((cHaveTcal = cPKSMS.sysCal().tableDesc().isColumn("TCAL"))) {
     167  if (cHaveSysCal && (cHaveTcal = cPKSMS.sysCal().tableDesc().isColumn("TCAL"))) {
    164168    cTcalCol.attach(cPKSMS.sysCal(), "TCAL");
    165169  }
    … …  
    176180  // Spectral data should always be present.
    177181  haveSpectra = True;
    178   cFloatDataCol.reference(msCols.floatData());
     182  cHaveDataCol = False;
     183  cHaveCorrectedDataCol = False;
     184  //String telName = antennaCols.name()(0);
     185  ROMSObservationColumns observationCols(cPKSMS.observation());
     186  String telName = observationCols.telescopeName()(0);
     187  //cATF = (telName.contains("DA41") || telName.contains("DV01"));
     188  cATF = telName.contains("ATF");
     189
     190  if (cHaveDataCol = cPKSMS.isColumn(MSMainEnums::DATA)) {
     191    if (cATF) {
     192      //try to read a single baseline interferometeric data
     193      //and treat it as single dish data
     194      //maybe extended for ALMA commissioning later
     195      cDataCol.reference(msCols.data());
     196      if (cHaveCorrectedDataCol = cPKSMS.isColumn(MSMainEnums::CORRECTED_DATA)) {
     197        //cerr<<"Do have CORRECTED_DATA column"<<endl;
     198        cCorrectedDataCol.reference(msCols.correctedData());
     199      }
     200    }
     201  }
     202  else {
     203    cFloatDataCol.reference(msCols.floatData());
     204  }
    179205  cFlagCol.reference(msCols.flag());
    180206
    181207
    182   if ((cGetXPol = cPKSMS.isColumn(MSMainEnums::DATA))) {
     208  if (cGetXPol = (cPKSMS.isColumn(MSMainEnums::DATA) && (!cATF))) {
    183209    if ((cHaveXCalFctr = cPKSMS.tableDesc().isColumn("XCALFCTR"))) {
    184210      cXCalFctrCol.attach(cPKSMS, "XCALFCTR");
    … …  
    223249  haveXPol = False;
    224250
    225   if (cGetXPol) {
     251  if (cGetXPol && !(cATF)) {
    226252    for (Int irow = 0; irow < cNRow; irow++) {
    227253      if (cDataCol.isDefined(irow)) {
    … …  
    303329
    304330  fluxunit = "";
    305   const TableRecord& keywordSet
    306      = cFloatDataCol.columnDesc().keywordSet();
     331  if (cHaveDataCol) {
     332    const TableRecord& keywordSet2
     333       = cDataCol.columnDesc().keywordSet();
     334    if(keywordSet2.isDefined("UNIT")) {
     335      fluxunit = keywordSet2.asString("UNIT");
     336    }
     337  } else {
     338    const TableRecord& keywordSet
     339       = cFloatDataCol.columnDesc().keywordSet();
     340    if(keywordSet.isDefined("UNIT")) {
     341      fluxunit = keywordSet.asString("UNIT");
     342    }
     343  }
     344
     345/***
     346  const TableRecord& keywordSet
     347       = cFloatDataCol.columnDesc().keywordSet();
    307348  if(keywordSet.isDefined("UNIT")) {
    308349    fluxunit = keywordSet.asString("UNIT");
    309350  }
    310 
     351***/
    311352  // Coordinate equinox.
    312353  ROMSPointingColumns pointingCols(cPKSMS.pointing());
    … …  
    529570        Double          &bandwidth,
    530571        Double          &freqInc,
    531         Double          &restFreq,
     572        Vector<Double>  &restFreq,
    532573        Vector<Float>   &tcal,
    533574        String          &tcalTime,
    … …  
    584625    }
    585626  }
    586 
    587627  // Renumerate scan no. Here still is 1-based
    588628  //scanNo = cScanNoCol(cIdx) - cScanNoCol(0) + 1;
    … …  
    624664
    625665  // Systemic velocity.
    626   if (!cHaveSrcVel) {
     666  if (!cHaveSrcVel || cATF) {
    627667    srcVel = 0.0f;
    628668  } else {
    … …  
    634674  Bool cGBT = telescope.contains("GBT");
    635675  // Observation type.
    636   Int stateId = cStateIdCol(cIdx);
    637   if (stateId == -1) {
    638     //obsMode = " ";
    639   } else {
    640     obsMode = cObsModeCol(stateId);
    641     Bool sigState =cSigStateCol(stateId);
    642     Bool refState =cRefStateCol(stateId);
    643     //DEBUG
    644     //cerr <<"stateid="<<stateId<<" obsmode="<<obsMode<<endl;
    645     if (cGBT) {
    646       // split the obsMode string and append a proper label
    647       // (these are GBT specific)
    648       int epos = obsMode.find_first_of(':');
    649       int nextpos = obsMode.find_first_of(':',epos+1);
    650       string obsMode1 = obsMode.substr(0,epos);
    651       string obsMode2 = obsMode.substr(epos+1,nextpos-epos-1);
     676  // check if State Table exist
     677  //Bool cHaveStateTab=Table::isReadable(cPKSMS.stateTableName());
     678  Int stateId = 0;
     679  Int StateNRow = 0;
     680  StateNRow=cObsModeCol.nrow();
     681  if (Table::isReadable(cPKSMS.stateTableName())) {
     682        obsMode = " ";
     683    if (StateNRow > 0) {
     684      stateId = cStateIdCol(cIdx);
     685      if (stateId == -1) {
     686        //obsMode = " ";
     687      } else {
     688        obsMode = cObsModeCol(stateId);
     689        Bool sigState =cSigStateCol(stateId);
     690        Bool refState =cRefStateCol(stateId);
     691        //DEBUG
     692        //cerr <<"stateid="<<stateId<<" obsmode="<<obsMode<<endl;
     693        if (cGBT) {
     694          // split the obsMode string and append a proper label
     695          // (these are GBT specific)
     696          int epos = obsMode.find_first_of(':');
     697          int nextpos = obsMode.find_first_of(':',epos+1);
     698          string obsMode1 = obsMode.substr(0,epos);
     699          string obsMode2 = obsMode.substr(epos+1,nextpos-epos-1);
    652700     
    653       //cerr <<"obsMode2= "<<obsMode2<<endl;
    654       if (!srcName.contains("_ps")
    655           &&!srcName.contains("_psr")
    656           &&!srcName.contains("_nod")
    657           &&!srcName.contains("_fs")
    658           &&!srcName.contains("_fsr")) {
    659         // if Nod mode observation , append '_nod'
    660         if (obsMode1 == "Nod") {
    661           srcName.append("_nod");
    662         } else if (obsMode1 == "OffOn") {
    663         // for GBT position switch observations (OffOn or OnOff)
    664           if (obsMode2 == "PSWITCHON") srcName.append("_ps");
    665           if (obsMode2 == "PSWITCHOFF") srcName.append("_psr");
    666         } else {
    667           if (obsMode2 == "FSWITCH") {
    668           // for GBT frequency switch mode
    669             if (sigState) srcName.append("_fs");
    670             if (refState) srcName.append("_fsr");
     701          //cerr <<"obsMode2= "<<obsMode2<<endl;
     702          if (!srcName.contains("_ps")
     703              &&!srcName.contains("_psr")
     704              &&!srcName.contains("_nod")
     705              &&!srcName.contains("_fs")
     706              &&!srcName.contains("_fsr")) {
     707            // if Nod mode observation , append '_nod'
     708            if (obsMode1 == "Nod") {
     709              srcName.append("_nod");
     710            } else if (obsMode1 == "OffOn") {
     711            // for GBT position switch observations (OffOn or OnOff)
     712              if (obsMode2 == "PSWITCHON") srcName.append("_ps");
     713              if (obsMode2 == "PSWITCHOFF") srcName.append("_psr");
     714            } else {
     715              if (obsMode2 == "FSWITCH") {
     716              // for GBT frequency switch mode
     717                if (sigState) srcName.append("_fs");
     718                if (refState) srcName.append("_fsr");
     719              }
     720            }
    671721          }
    672         } 
    673       }
    674     }
    675   } 
     722        }
     723      } 
     724    }
     725  }
    676726  // CAL state
    677727  // this should be apply just for GBT data?
    678728  Double Cal;
    679   if (stateId==-1) {
     729  if (stateId==-1 || StateNRow==0) {
    680730    Cal = 0;
    681731  } else {
    … …  
    694744  Vector<Double> chanFreq = cChanFreqCol(iIF);
    695745  if (nChan == 1) {
    696     cout << "The input is continuum data. "<< endl;
    697746    freqInc  = chanFreq(0);
    698747    refFreq  = chanFreq(0);
    699748    restFreq = 0.0f;
    700749  } else {
     750 
    701751    if (cStartChan(iIF) <= cEndChan(iIF)) {
    702752      freqInc = chanFreq(1) - chanFreq(0);
    … …  
    705755    }
    706756    refFreq  = chanFreq(cRefChan(iIF)-1);
    707     restFreq = cSrcRestFrqCol(srcId)(IPosition(1,0));
     757    Bool HaveSrcRestFreq= cSrcRestFrqCol.isDefined(srcId);
     758    if (HaveSrcRestFreq) {
     759      //restFreq = cSrcRestFrqCol(srcId)(IPosition(1,0));
     760      restFreq = cSrcRestFrqCol(srcId);
     761    } else {
     762      restFreq = 0.0f;
     763    }
    708764  }
    709765  bandwidth = abs(freqInc * nChan);
    … …  
    720776
    721777  // Find the appropriate entry in the WEATHER subtable.
    722   Vector<Double> wTimes = cWeatherTimeCol.getColumn();
    723   Int weatherIdx;
    724   for (weatherIdx = wTimes.nelements()-1; weatherIdx >= 0; weatherIdx--) {
    725     if (cWeatherTimeCol(weatherIdx) <= time) {
    726       break;
    727     }
    728   }
    729 
    730   if (weatherIdx < 0) {
     778  //Bool cHaveStateTab=Table::isReadable(cPKSMS.stateTableName());
     779  Bool cHaveWeatherTab = Table::isReadable(cPKSMS.weatherTableName());
     780  Int weatherIdx=-1;
     781  if (cHaveWeatherTab) {
     782    Vector<Double> wTimes = cWeatherTimeCol.getColumn();
     783    for (weatherIdx = wTimes.nelements()-1; weatherIdx >= 0; weatherIdx--) {
     784      if (cWeatherTimeCol(weatherIdx) <= time) {
     785        break;
     786      }
     787    }
     788  }
     789  if (weatherIdx < 0 || !cHaveWeatherTab) {
    731790    // No appropriate WEATHER entry.
    732791    pressure    = 0.0f;
    … …  
    874933    baseSub.resize(0,0);
    875934  }
    876 
    877935  // Get spectral data.
    878936  if (cGetSpectra) {
    879937    Matrix<Float> tmpData;
    880938    Matrix<Bool>  tmpFlag;
    881     cFloatDataCol.getSlice(cIdx, cDataSel(iIF), tmpData, True);
     939    if (cHaveDataCol) {
     940      Matrix<Complex> tmpCmplxData;
     941      Matrix<Float> tmpReData;
     942      Matrix<Float> tmpImData;
     943      //cerr<<"reading spectra..."<<endl;
     944      //# TODO - should have a flag to user to select DATA or CORRECTED_DATA
     945      //# currently just automatically determined, --- read CORRECTED one
     946      //# if the column exist.
     947      if (cHaveCorrectedDataCol) {
     948        cCorrectedDataCol.getSlice(cIdx, cDataSel(iIF), tmpCmplxData, True);
     949      } else {
     950        cDataCol.getSlice(cIdx, cDataSel(iIF), tmpCmplxData, True);
     951      }
     952      tmpReData = real(tmpCmplxData);
     953      tmpImData = imag(tmpCmplxData);
     954      tmpData = sqrt(tmpReData*tmpReData + tmpImData*tmpImData);
     955    } else {
     956      cFloatDataCol.getSlice(cIdx, cDataSel(iIF), tmpData, True);
     957    }
    882958    cFlagCol.getSlice(cIdx, cDataSel(iIF), tmpFlag, True);
    883959
    … …  
    909985  // Get cross-polarization data.
    910986  if (cGetXPol) {
     987    //cerr<<"cGetXPol="<<cGetXPol<<endl;
     988    //cerr<<"cHaveXCalFctr="<<cHaveXCalFctr<<endl;
     989
    911990    if (cHaveXCalFctr) {
    912991      cXCalFctrCol.get(cIdx, xCalFctr);
    … …  
    915994    }
    916995
    917     cDataCol.get(cIdx, xPol, True);
    918 
    919     if (cEndChan(iIF) < cStartChan(iIF)) {
    920       Complex ctmp;
    921       Int jchan = nChan - 1;
    922       for (Int ichan = 0; ichan < nChan/2; ichan++, jchan--) {
    923         ctmp = xPol(ichan);
    924         xPol(ichan) = xPol(jchan);
    925         xPol(jchan) = ctmp;
    926       }
    927     }
    928   }
    929 
     996    if(!cATF) {
     997      cDataCol.get(cIdx, xPol, True);
     998
     999      if (cEndChan(iIF) < cStartChan(iIF)) {
     1000        Complex ctmp;
     1001        Int jchan = nChan - 1;
     1002        for (Int ichan = 0; ichan < nChan/2; ichan++, jchan--) {
     1003          ctmp = xPol(ichan);
     1004          xPol(ichan) = xPol(jchan);
     1005          xPol(jchan) = ctmp;
     1006        }
     1007      }
     1008    }
     1009  }
     1010  /**
     1011  cerr<<"scanNo="<<scanNo<<endl;
     1012  cerr<<"cycleNo="<<cycleNo<<endl;
     1013  cerr<<"mjd="<<mjd<<endl;
     1014  cerr<<"interval="<<interval<<endl;
     1015  cerr<<"fieldName="<<fieldName<<endl;
     1016  cerr<<"srcNmae="<<srcName<<endl;
     1017  cerr<<"srcDir="<<srcDir<<endl;
     1018  cerr<<"srcPM="<<srcPM<<endl;
     1019  cerr<<"srcVel="<<srcVel<<endl;
     1020  cerr<<"obsMode="<<obsMode<<endl;
     1021  cerr<<"IFno="<<IFno<<endl;
     1022  cerr<<"refFreq="<<refFreq<<endl;
     1023  cerr<<"tcal="<<tcal<<endl;
     1024  cerr<<"direction="<<direction<<endl;
     1025  cerr<<"scanRate="<<scanRate<<endl;
     1026  cerr<<"tsys="<<tsys<<endl;
     1027  cerr<<"sigma="<<sigma<<endl;
     1028  cerr<<"calFctr="<<calFctr<<endl;
     1029  cerr<<"baseLin="<<baseLin<<endl;
     1030  cerr<<"baseSub="<<baseSub<<endl;
     1031  cerr<<"spectra="<<spectra.shape()<<endl;
     1032  cerr<<"flagged="<<flagged.shape()<<endl;
     1033  cerr<<"xCalFctr="<<xCalFctr<<endl;
     1034  cerr<<"xPol="<<xPol<<endl;
     1035  **/
    9301036  cIdx++;
    9311037
    … …  
    10031109    Matrix<Float> tmpData;
    10041110    Matrix<Bool>  tmpFlag;
    1005     cFloatDataCol.getSlice(cIdx, cDataSel(iIF), tmpData, True);
     1111    if (cHaveDataCol) {
     1112      Matrix<Complex> tmpCmplxData;
     1113      cDataCol.getSlice(cIdx, cDataSel(iIF), tmpCmplxData, True);
     1114      tmpData = real(tmpCmplxData);
     1115    } else {
     1116      cFloatDataCol.getSlice(cIdx, cDataSel(iIF), tmpData, True);
     1117    }
    10061118    cFlagCol.getSlice(cIdx, cDataSel(iIF), tmpFlag, True);
    10071119
Note: See TracChangeset for help on using the changeset viewer.