Changeset 16 for trunk/src


Ignore:
Timestamp:
06/30/04 11:00:23 (20 years ago)
Author:
mmarquar
Message:

Updated data container. Changed the axis order in the spectrum/flag arrays to [nBeam,nIF,nPol,nChan]

Location:
trunk/src
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • trunk/src/SDContainer.cc

    r14 r16  
    3838using namespace atnf_sd;
    3939
    40 SDContainer::SDContainer(uInt nBeam, uInt nIF, uInt nChan, uInt nPol)
     40SDContainer::SDContainer(uInt nBeam, uInt nIF, uInt nPol, uInt nChan)
    4141  : nBeam_(nBeam),
    42   nIF_(nIF),
    43   nChan_(nChan),
    44   nPol_(nPol),
    45   spectrum_(IPosition(4,nBeam,nIF,nChan,nPol)),
    46   flags_(IPosition(4,nBeam,nIF,nChan,nPol)),
    47   tsys_(IPosition(4,nBeam,nIF,nChan,nPol)) {
     42    nIF_(nIF),
     43    nPol_(nPol),
     44    nChan_(nChan),
     45    spectrum_(IPosition(4,nBeam,nIF,nPol,nChan)),
     46    flags_(IPosition(4,nBeam,nIF,nPol,nChan)),
     47    tsys_(IPosition(4,nBeam,nIF,nPol,nChan)) {
     48  uChar x = 0;
     49  flags_ = ~x;
     50}
     51
     52SDContainer::SDContainer(IPosition shp)
     53  : nBeam_(shp(0)),
     54    nIF_(shp(1)),
     55    nPol_(shp(2)),
     56    nChan_(shp(3)),
     57    spectrum_(shp),
     58    flags_(shp),
     59    tsys_(shp) {
    4860  uChar x = 0;
    4961  flags_ = ~x;
    … …  
    7284 
    7385  //Vector<Float> pols(nPol);
    74   ArrayAccessor<Float, Axis<0> > j(spec);
     86  ArrayAccessor<Float, Axis<1> > j(spec);
    7587  IPosition shp0 = spectrum_.shape();
    7688  IPosition shp1 = spec.shape();
    77   if ( (shp0(2) != shp1(0)) || (shp0(3) != shp1(1)) ) {
     89  if ( (shp0(2) != shp1(1)) || (shp0(3) != shp1(0)) ) {
    7890    cerr << "Arrays not conformant" << endl;     
    7991    return False;
    8092  }
    81   //uInt chani = 0;
    82   //uInt poli = 0;
    8393  // assert dimensions are the same....
    8494  for (ArrayAccessor<Float, Axis<2> > i(aa1);i != i.end(); ++i) {
    85     //pols = spec.row(chani);   
    86     ArrayAccessor<Float, Axis<1> > jj(j);
     95    ArrayAccessor<Float, Axis<0> > jj(j);
    8796    for (ArrayAccessor<Float, Axis<3> > ii(i);ii != ii.end(); ++ii) {
    88       //(*ii) = pols[poli];     
    8997      (*ii) = (*jj);
    90       //poli++;
    9198      jj++;
    9299    }
    93     //poli = 0;
    94     //chani++;
    95100    j++;
    96101  }
    … …  
    113118  aa1.reset(aa1.begin(whichIF));
    114119 
    115   ArrayAccessor<uChar, Axis<0> > j(flag);
     120  ArrayAccessor<uChar, Axis<1> > j(flag);
    116121  IPosition shp0 = flags_.shape();
    117122  IPosition shp1 = flag.shape();
    118   if ( (shp0(2) != shp1(0)) || (shp0(3) != shp1(1)) ) {
     123  if ( (shp0(2) != shp1(1)) || (shp0(3) != shp1(0)) ) {
    119124    cerr << "Arrays not conformant" << endl;     
    120125    return False;
    … …  
    123128  // assert dimensions are the same....
    124129  for (ArrayAccessor<uChar, Axis<2> > i(aa1);i != i.end(); ++i) {
    125     ArrayAccessor<uChar, Axis<1> > jj(j);
     130    ArrayAccessor<uChar, Axis<0> > jj(j);
    126131    for (ArrayAccessor<uChar, Axis<3> > ii(i);ii != ii.end(); ++ii) {
    127132      (*ii) = (*jj);
    … …  
    130135    j++;
    131136  }
     137  return True;
    132138}
    133139
    … …  
    139145  aa1.reset(aa1.begin(whichIF));
    140146  // assert dimensions are the same....
    141   uInt idx = 0;
    142 
    143 
    144   for (ArrayAccessor<Float, Axis<2> > i(aa1);i != i.end(); ++i) {   
    145     idx = 0;
     147  for (ArrayAccessor<Float, Axis<3> > i(aa1);i != i.end(); ++i) {   
    146148    ArrayAccessor<Float, Axis<0> > j(tsys);
    147     for (ArrayAccessor<Float, Axis<3> > ii(i);ii != ii.end(); ++ii) {
     149    for (ArrayAccessor<Float, Axis<2> > ii(i);ii != ii.end(); ++ii) {
    148150      (*ii) = (*j);
    149151      j++;

  • trunk/src/SDContainer.h

    r7 r16  
    6464
    6565public:
    66   SDContainer(uInt nBeam, uInt nIF, uInt nChan, uInt nPol);
     66  SDContainer(uInt nBeam, uInt nIF, uInt nPol, uInt nChan);
     67  SDContainer(IPosition shp);
    6768
    6869  virtual ~SDContainer();
    … …  
    9495 
    9596private:
    96   uInt nBeam_,nIF_,nChan_,nPol_;
     97  uInt nBeam_,nIF_,nPol_,nChan_;
    9798
    98   // (nBeam,nIF,nChannel,nPol)
     99  // (nBeam,nIF,nPol,nChannel)
    99100  Array<Float>    spectrum_; 
    100101  Array<uChar>    flags_;
    101   // (nBeam,nIF,[nChannel],nPol) Tsys is not really a function of
     102  // (nBeam,nIF,nPol,[nChannel]) Tsys is not really a function of
    102103  // channel, but this makes it easier to work with at the expense of
    103104  // a little memory

  • trunk/src/SDMemTable.cc

    r8 r16  
    3030//#---------------------------------------------------------------------------
    3131
     32
    3233#include <aips/iostream.h>
    3334#include <aips/Arrays/Array.h>
    … …  
    6061}
    6162
    62 SDMemTable::SDMemTable(const SDMemTable& other) {
     63SDMemTable::SDMemTable(const SDMemTable& other, Bool clear) {
    6364  this->IFSel_= other.IFSel_;
    6465  this->beamSel_= other.beamSel_;
    … …  
    6869  this->table_ = other.table_.copyToMemoryTable(String("dummy"));
    6970  // clear all rows()
    70   this->table_.removeRow(this->table_.rowNumbers());
     71  if (clear) {
     72    this->table_.removeRow(this->table_.rowNumbers());
     73  } else {
     74    this->IFSel_ = other.IFSel_;
     75    this->beamSel_ = other.beamSel_;
     76    this->polSel_ = other.polSel_;
     77  }
    7178}
    7279
    … …  
    101108  td.addColumn(ArrayColumnDesc<Float>("TSYS")); 
    102109  td.addColumn(ScalarColumnDesc<Int>("SCANID")); 
     110  td.addColumn(ScalarColumnDesc<Double>("INTERVAL")); 
    103111  // Now create a new table from the description.
    104112  SetupNewTable aNewTab(name_, td, Table::New);
    … …  
    147155}
    148156
    149 bool SDMemTable::setChannels(const std::vector<int>& whichChans) {
    150   //std::vector<bool>::iterator it;
    151   //std::fill(chanMask_.begin(), chanMask_.end(), false);
    152   //for (it = whichChans.begin(); it != whichChans.end(); it++) {
    153   //chanMask_[*it] = true;
    154   //}
    155   cout << "SDMemTable::setChannels() disabled" << endl;
     157bool SDMemTable::setMask(const std::vector<int>& whichChans) {
     158  ROArrayColumn<uChar> spec(table_, "FLAGTRA");
     159 
     160  std::vector<int>::iterator it;
     161  uInt n = spec.shape(0)(3);
     162  chanMask_.resize(n,true);
     163  for (it = whichChans.begin(); it != whichChans.end(); it++) {
     164    if (*it < n)
     165      chanMask_[*it] = false;
     166  }
    156167  return true;
    157168}
    158169
    159 std::vector<bool> SDMemTable::getMask() const {
    160   return chanMask_;
     170std::vector<bool> SDMemTable::getMask(Int whichRow) const {
     171  std::vector<bool> mask;
     172  ROArrayColumn<uChar> spec(table_, "FLAGTRA");
     173  Array<uChar> arr;
     174  spec.get(whichRow, arr);
     175  ArrayAccessor<uChar, Axis<0> > aa0(arr);
     176  aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
     177  ArrayAccessor<uChar, Axis<1> > aa1(aa0);
     178  aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
     179  ArrayAccessor<uChar, Axis<2> > aa2(aa1);
     180  aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
     181
     182  Bool useUserMask = ( chanMask_.size() == arr.shape()(3) );
     183
     184  std::vector<bool> tmp;
     185  tmp = chanMask_; // WHY the fxxx do I have to make a copy here
     186  std::vector<bool>::iterator miter;
     187  miter = tmp.begin();
     188
     189  for (ArrayAccessor<uChar, Axis<3> > i(aa2); i != i.end(); ++i) {
     190    bool out =!static_cast<bool>(*i);
     191    if (useUserMask) {
     192      out = out && (*miter);
     193      miter++;
     194    }
     195    mask.push_back(out);
     196  } 
     197  return mask;
    161198}
    162199std::vector<float> SDMemTable::getSpectrum(Int whichRow) const {
    … …  
    170207  ArrayAccessor<Float, Axis<1> > aa1(aa0);
    171208  aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
    172   for (ArrayAccessor<Float, Axis<2> > i(aa1); i != i.end(); ++i) {
    173     ArrayAccessor<Float, Axis<3> > aa2(i);
    174     aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
    175     spectrum.push_back(*aa2);
     209  ArrayAccessor<Float, Axis<2> > aa2(aa1);
     210  aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
     211  for (ArrayAccessor<Float, Axis<3> > i(aa2); i != i.end(); ++i) {
     212    spectrum.push_back(*i);
    176213  }
    177214  return spectrum;
    178215}
    179216
    180 MaskedArray<Float> SDMemTable::rowAsMaskedArray(uInt whichRow) {
     217MaskedArray<Float> SDMemTable::rowAsMaskedArray(uInt whichRow,
     218                                                Bool useSelection) {
    181219  ROArrayColumn<Float> spec(table_, "SPECTRA");
    182220  Array<Float> arr;
    … …  
    185223  spec.get(whichRow, arr);
    186224  flag.get(whichRow, farr);
    187   IPosition test(4,0,0,100,0);
    188225  Array<Bool> barr(farr.shape());convertArray(barr, farr);
    189226  MaskedArray<Float> marr;
    190   marr.setData(arr,!barr);
     227  if (useSelection) {
     228    ArrayAccessor<Float, Axis<0> > aa0(arr);
     229    aa0.reset(aa0.begin(uInt(beamSel_)));//go to beam
     230    ArrayAccessor<Float, Axis<1> > aa1(aa0);
     231    aa1.reset(aa1.begin(uInt(IFSel_)));// go to IF
     232    ArrayAccessor<Float, Axis<2> > aa2(aa1);
     233    aa2.reset(aa2.begin(uInt(polSel_)));// go to pol
     234   
     235    ArrayAccessor<Bool, Axis<0> > baa0(barr);
     236    baa0.reset(baa0.begin(uInt(beamSel_)));//go to beam
     237    ArrayAccessor<Bool, Axis<1> > baa1(baa0);
     238    baa1.reset(baa1.begin(uInt(IFSel_)));// go to IF
     239    ArrayAccessor<Bool, Axis<2> > baa2(baa1);
     240    baa2.reset(baa2.begin(uInt(polSel_)));// go to pol
     241
     242    Vector<Float> a(arr.shape()(3));
     243    Vector<Bool> b(barr.shape()(3));
     244    ArrayAccessor<Float, Axis<0> > a0(a);
     245    ArrayAccessor<Bool, Axis<0> > b0(b);
     246
     247    ArrayAccessor<Bool, Axis<3> > j(baa2);
     248    for (ArrayAccessor<Float, Axis<3> > i(aa2); i != i.end(); ++i) {
     249      (*a0) = (*i);
     250      (*b0) = !(*j);
     251      j++;
     252      a0++;
     253      b0++;
     254    }
     255    marr.setData(a,b);
     256  } else {
     257    marr.setData(arr,!barr);
     258  }
    191259  return marr;
    192260}
    193 
    194261
    195262Float SDMemTable::getTsys(Int whichRow) const {
    … …  
    199266  Float out;
    200267  IPosition ip(arr.shape());
    201   ip(0) = beamSel_;ip(1) = IFSel_;ip(2)=0;ip(3) = polSel_;
     268  ip(0) = beamSel_;ip(1) = IFSel_;ip(2) = polSel_;ip(3)=0;
    202269  out = arr(ip);
    203270  return out;
    … …  
    211278  ArrayColumn<Float> ts(table_, "TSYS");
    212279  ScalarColumn<Int> scan(table_, "SCANID");
     280  ScalarColumn<Double> integr(table_, "INTERVAL");
    213281
    214282  uInt rno = table_.nrow();
    … …  
    221289  ts.put(rno, sdc.getTsys());
    222290  scan.put(rno, sdc.scanid);
     291  integr.put(rno, sdc.interval);
     292 
    223293  return true;
    224294}
    … …  
    241311      previous = current;     
    242312      count++;
    243       cout << count << "\t" << name << endl;
     313      cout << count << "\t"
     314           << name
     315           << endl;
    244316    }
    245317  }
    … …  
    247319  cout << "in " << count << "scans." << endl;
    248320}
    249 
     321/*
     322void maskChannels(const std::vector<Int>& whichChans ) {
     323 
     324  std::vector<int>::iterator it;
     325  ArrayAccessor<uChar, Axis<2> > j(flags_);
     326  for (it = whichChans.begin(); it != whichChans.end(); it++) {
     327    j.reset(j.begin(uInt(*it)));
     328    for (ArrayAccessor<uChar, Axis<0> > i(j); i != i.end(); ++i) {
     329      for (ArrayAccessor<uChar, Axis<1> > ii(i); ii != ii.end(); ++ii) {
     330        for (ArrayAccessor<uChar, Axis<3> > iii(ii);
     331             iii != iii.end(); ++iii) {   
     332          (*iii) =
     333        }
     334      }
     335    }
     336  }
     337 
     338}
     339*/

  • trunk/src/SDMemTable.h

    r2 r16  
    5151public:
    5252  SDMemTable(const std::string& name= "SDInputTable.tbl");
    53   SDMemTable(const SDMemTable& other);
     53  SDMemTable(const SDMemTable& other, Bool clear=False);
    5454
    5555  SDMemTable(const Table& tab, Int scanID);
    … …  
    6060 
    6161  virtual std::vector<float> getSpectrum(Int whichRow) const;
    62   virtual std::vector<bool> getMask() const;
     62  virtual std::vector<bool> getMask(Int whichRow) const;
    6363 
    64   MaskedArray<Float> rowAsMaskedArray(uInt whichRow);
     64  MaskedArray<Float> rowAsMaskedArray(uInt whichRow,
     65                                      Bool useSelection = False);
    6566
    6667  virtual Float getTsys(Int whichRow) const;
    … …  
    7778
    7879  //sets the mask
    79   virtual bool setChannels(const std::vector<int>& whichChans);
     80  virtual bool setMask(const std::vector<int>& whichChans);
    8081 
    8182  virtual void summary() const;

  • trunk/src/SDMemTableWrapper.h

    r10 r16  
    6161  double getTime(int whichRow) {return table_->getTime(whichRow);}
    6262
    63   std::vector<bool> getMask() const { return table_->getMask(); }
     63  std::vector<bool> getMask(int whichRow) const {
     64    return table_->getMask(whichRow);
     65  }
     66  bool setMask(const std::vector<int> mvals) const {
     67    return table_->setMask(mvals);
     68  }
    6469
    6570  std::string getSourceName(int whichRow) {
    … …  
    111116                                             off.getCP()));
    112117  }
     118  SDMemTableWrapper multiply(const SDMemTableWrapper& in,
     119                             Float factor) {
     120    return SDMemTableWrapper(SDMath::multiply(in.getCP(),factor));
     121  }
    113122};
    114123

  • trunk/src/SDReader.cc

    r14 r16  
    3030//#---------------------------------------------------------------------------
    3131#include <atnf/PKSIO/PKSreader.h>
    32 
     32#include <aips/Quanta/MVTime.h>
    3333#include "SDContainer.h"
    3434#include "SDReader.h"
    3535
    3636using namespace atnf_sd;
     37
     38void SDHeader::print() const {
     39  MVTime mvt(this->utc);
     40
     41  cout << "Observer: " << this->observer << endl
     42       << "Project: " << this->project << endl
     43       << "Obstype: " << this->obstype << endl
     44       << "Antenna: " << this->antennaname << endl
     45       << "Ant. Position: " << this->antennaposition << endl
     46       << "Equinox: " << this->equinox << endl
     47       << "Freq. ref.: " << this->freqref << endl
     48       << "Ref. frequency: " << this->reffreq << endl
     49       << "Bandwidth: "  << this->bandwidth << endl
     50       << "Time (utc): "
     51       << mvt.string()
     52       << endl;
     53  //setprecision(10) << this->utc << endl;
     54}
     55
    3756
    3857SDReader::SDReader() :
    … …  
    85104 
    86105  // Get basic parameters.
    87   Double bandwidth, refFreq;
    88106  header_ = SDHeader();
    89107  header_.nchan = nChan_;
    … …  
    152170  while ( (cursor_ <= seq[n-1]) || getAll) {
    153171    // only needs to be create if in seq
    154     SDContainer sc(header_.nbeam,header_.nif,header_.nchan,header_.npol);
     172    SDContainer sc(header_.nbeam,header_.nif,header_.npol,header_.nchan);
    155173
    156174    // iterate over one correlator integration cycle = nBeam*nIF
    157175    for (uInt row=0; row < stepsize; row++) {
    158 
    159       //cerr << "SDReader starting reading process row = " << row << endl;
    160 
    161176      // add scanid from GROUP field -- this will remove the need for
    162177      // stepsize as well
    … …  
    171186                             tsys, sigma, calFctr, baseLin, baseSub,
    172187                             spectra, flagtra, xCalFctr, xPol);         
    173 
    174       // cerr << "SDReader row  read"<< endl;
    175 
    176188      if (status) {
    177189        if (status == -1) {
    … …  
    183195      }     
    184196      // if in the given list
    185       if ( cursor_ == seq[seqi]) {
    186 
    187         //cerr << "SDReader cursor == seq[i]" << endl;
    188        
     197      if (cursor_ == seq[seqi]) {
    189198        // add integration cycle
    190199        if (row==0) {
    … …  
    193202          //focusRot, temperature, pressure, humidity, windSpeed,
    194203          //windAz  srcDir, srcPM, srcVel
     204          sc.timestamp = mjd;
     205          sc.interval = interval;         
     206          sc.sourcename = srcName;
    195207        }
    196 
    197         // cerr << "SDReader::read -  before SDContainer" << endl;
    198 
    199208        // add specific info
    200209        // IFno beamNo are 1-relative
    201         // refPix = nChan/2+1 in  Integer arith.!
    202        
    203         uInt refPix = header_.nchan/2+1;
     210        // refPix = nChan/2+1 in  Integer arith.!       
     211        //uInt refPix = header_.nchan/2+1;
    204212        //uInt frqslot = sdft.addFrequency(refPix, refFreq, freqInc);
    205213
    … …  
    208216          prevName = srcName;//temp
    209217        }//temp
    210 
    211         sc.sourcename = srcName;
    212         sc.timestamp = mjd;
    213         sc.interval = interval;
    214218        sc.scanid = scanid;
    215219        //sc.setFrequencyMap(frqslot,IFno-1);
    … …  
    218222        sc.setTsys(tsys, beamNo-1, IFno-1);
    219223        //sc.addPointing(direction, beamNo-1);
    220 
    221         // cerr << "SDReader::read -  after SDContainer" << endl;
    222        
    223224      }
    224225    }
    … …  
    226227      // insert container into table/list
    227228      table_->putSDContainer(sc);
    228       //cout << "Reading integration " << seqi << endl;
    229229      seqi++;// next in list
    230230    }
    231     cursor_++;
     231    cursor_++;// increment position in file
    232232  }
    233233  return status;

  • trunk/src/SDReader.h

    r2 r16  
    6363  Double bandwidth;
    6464  Double utc;
    65   void print() {
    66     cout << "Observer: " << this->observer << endl
    67          << "Project: " << this->project << endl
    68          << "Obstype: " << this->obstype << endl
    69          << "Antenna: " << this->antennaname << endl
    70          << "Ant. Position: " << this->antennaposition << endl
    71          << "Equinox: " << this->equinox << endl
    72          << "Freq. ref.: " << this->freqref << endl
    73          << "Ref. frequency: " << this->reffreq << endl
    74          << "Bandwidth: "  << this->bandwidth << endl
    75          << "Time (utc): " << setprecision(10) << this->utc << endl;
    76   };
     65  void print() const ;
    7766};
    7867
    … …  
    10089
    10190private:
    102   Int nBeam_,nIF_,nChan_,nPol_;
     91  Int nBeam_,nIF_,nPol_,nChan_;
    10392  Bool getHeader();
    10493  PKSreader* reader_; 

  • trunk/src/SDTemplates.cc

    r11 r16  
    5555template LogicalArray operator>=(Array<Float> const &, Float const &);
    5656template Array<Float>& operator/=(Array<Float>&, MaskedArray<Float> const&);
    57 //template Array<Float> const& operator*=(Array<Float> const&, MaskedArray<Float> const&);
     57template MaskedArray<Float> const& operator*=(MaskedArray<Float> const&, Float const&);
    5858template MaskedArray<Float> operator-(MaskedArray<Float> const&, MaskedArray<Float> const&);

  • trunk/src/python_SDMemTable.cc

    r2 r16  
    3232
    3333#include <boost/python.hpp>
    34 
     34#include <boost/python/args.hpp>
    3535#include "SDMemTableWrapper.h"
    3636
    … …  
    4646    .def("getscan", &SDMemTableWrapper::getScan)
    4747    .def("getspectrum", &SDMemTableWrapper::getSpectrum)
    48     //.def("getmask", &SDMemTableWrapper::getMask)
     48    .def("getmask", &SDMemTableWrapper::getMask)
    4949    .def("gettsys", &SDMemTableWrapper::getTsys)
    5050    .def("getsourcename", &SDMemTableWrapper::getSourceName)
    … …  
    5353    .def("getbeam", &SDMemTableWrapper::getBeam)
    5454    .def("getpol", &SDMemTableWrapper::getPol)
    55     .def("setif", &SDMemTableWrapper::setIF)
     55    .def("setif", &SDMemTableWrapper::setIF,
     56         (boost::python::arg("whichIF")=0) )
    5657    .def("setbeam", &SDMemTableWrapper::setBeam)
    5758    .def("setpol", &SDMemTableWrapper::setPol)
     59    .def("setmask", &SDMemTableWrapper::setMask)
    5860    .def("makepersistent",  &SDMemTableWrapper::makePersistent)
    5961    .def("summary",  &SDMemTableWrapper::summary)
Note: See TracChangeset for help on using the changeset viewer.