Importmiriad.h

Go to the documentation of this file.

//# Importmiriad: miriad dataset to MeasurementSet conversion
//# Copyright (C) 1997,2000,2001,2002,2013,2015
//# Associated Universities, Inc. Washington DC, USA.
//#
//# This program is free software; you can redistribute it and/or modify
//# it under the terms of the GNU General Public License as published by
//# the Free Software Foundation; either version 2 of the License, or
//# (at your option) any later version.
//#
//# This program is distributed in the hope that it will be useful,
//# but WITHOUT ANY WARRANTY; without even the implied warranty of
//# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//# GNU General Public License for more details.
//# 
//# You should have received a copy of the GNU General Public License
//# along with this program; if not, write to the Free Software
//# Foundation, Inc., 675 Masve, Cambridge, MA 02139, USA.
//#

#ifndef MIRIAD_IMPORTMIRIAD_H
#define MIRIAD_IMPORTMIRIAD_H

//#Includes

#include <casa/Inputs/Input.h>
#include <casa/OS/File.h>
#include <casa/Utilities/GenSort.h>

#include <casa/Arrays/Cube.h>
#include <casa/Arrays/Matrix.h>
#include <casa/Arrays/Vector.h>
#include <casa/Arrays/ArrayMath.h>
#include <casa/Arrays/ArrayUtil.h>
#include <casa/Arrays/ArrayLogical.h>
#include <casa/Arrays/MatrixMath.h>

#include <measures/Measures.h>
#include <measures/Measures/MPosition.h>
#include <measures/Measures/MeasData.h>
#include <measures/Measures/Stokes.h>

#include <tables/Tables.h>
#include <tables/Tables/TableInfo.h>

#include <ms/MeasurementSets.h> 

#include <mirlib/maxdimc.h>
#include <mirlib/miriad.h>

#include <casa/namespace.h>
namespace casa { //# NAMESPACE CASA - BEGIN


// MIRIAD dataset MeasurementSet filler. Derived from carmafiller.
//
// Despite that this program is written in C++ and uses classes, it is
// really a set of routines manipulating a huge common block.....
//
// Acknowledgement: this program was originally cloned off uvfitsfiller.cc
//                  whose code is now in MSFitsInput.cc
//                  Ray Plante wrote the other filler, with which we had
//                  some cross talk
// 
//
// Limitations, Caveats and Remaining ToDo's
//    - does not apply the various (miriad) complex gain corrections
//      (uvcat and friends should be used for that in miriad, it would be
//      silly to duplicate that code here)
//    - this code won't run if Double!=double or Float!=float or Int!=int
//      On some future 64bit machines this may cause problems?
//    - CARMA type arrays when dishes are not the same size
//    - do not mix arrays from different telescopes (e.g. BIMA and ATCA)
//      but it should handle different configurations (ARRAYs) of the same
//      instrument
//    - handle multiple array configuration datasets that were UVCAT'ed
//      (needs to count antennae up at each array configuration)
//    - true miriad storage manager ?
//
//    - check UV override, is that being handled transparently?
//
//    - spectral windows layout  (code present, fix table access descriptors) 
//      CHECK: what happens if e.g. uvcat does some line= preprocessing
//    - correct restfreq's in source tables that relate to the spectral windows
//      (this code currently coredumps)
//    - make win=,narrow=,wide=0 cause non of them to be written in that selection
//    - add weather table
//    - although no more wides are written, there is a way to make a small MS
//      from the wides only by preprocessing in miriad, but conceivably
//      could be done here as well??  -- except the narrow= keyword seems to write
//      files that suffer from the Table array conformance error
//
//      Todo: to deal with multiple zoom setups (freq changes), we need
//      to track which spectral windows each source appears to make a proper
//      source table indexed by source id and spw id


//  History:
//   Spring 1997:   written (cloned off uvfitsfiller)          Peter Teuben
//   July 1997:     Y2K, fixed table-interface                          PJT
//   Dec 1997:      fixed wideband only data (e.g. uvgen)               PJT
//   ???            somebody fixed up this code for some new release    ???
//   May 2000:      fixed up for various new AIPS++ conventions         PJT
//                  and added multi-source & field                      PJT
//   Sep 2000:      development now on linux, converted to OldMS        PJT
//   Dec 2000:      Conversion to MS (MS2)                              PJT
//                  typical compile time: (P600/256M/15MBps HD: 36")
//                  typical 3c273 conversion time:  2.6" (5.3->11.3 MB)
//                  Cf. that to "uvio" processing time, which runs at disk I/O
//                  speed (15 MB/s on the beforementioned laptop)
//   Jan-Feb 2001:  window layout, restfreq's and various things to get 
//                  msmultiscale to work; added some syscal support     PJT
//   April 2001:    final cleanup for the 1.5 AIPS++ release -          PJT
//                  [note that from this moment on bimafiller is for
//                   experimental use and deprecated, the new version
//                   will be called 'mirfiller'
//   Oct 2001:      mirlib changed location
//   Sep 2009:      revived in CASA, added various tables for CARMA
//                  also renamed from bimafiller to carmafiller
//   Sep 2011:      initial release
//   Mar 2013:      Rename to importmiriad and make it work for ATCA 
//                  CABB data (more channels, windows and 4 pols)       MHW
//   May 2014:      Fix some compiler warnings, add TOPO option         MHW
//   Jul 2015:      Cope with data with changing freq setups,e.g.,
//                  CABB zoom data with multiple sources                MHW
//   Jan 2016       Turn from standalone app into tool+task             MHW

// a placeholder for the MIRIAD spectral window configuration
// see also the MIRIAD programmers guide, and its appendix on UV Variables 
// for CARMA we have MAXWIN = MAXWIDE (16 normally, 6 in 2006-2009, more later)
// narrow band (plus MAXWIDE space to point into the wide band data too)

typedef struct window {     // CASA defines everything mid-band mid-interval
  int    nspect;                   // number of valid windows (<=MAXWIN, typically 16)
  int    nschan[MAXWIN+MAXWIDE];   // number of channels in a window
  int    ischan[MAXWIN+MAXWIDE];   // starting channel of a window (1-based)
  double sdf[MAXWIN+MAXWIDE];      // channel separation 
  double sfreq[MAXWIN+MAXWIDE];    // frequency of first channel in window (doppler changes)
  double restfreq[MAXWIN+MAXWIDE]; // rest freq, if appropriate
  char   code[MAXWIN+MAXWIDE];     // code to CASA identification (N, W or S; S not used anymore)

  // wide band (for CARMA these are the spectral window averages - i.e. nspect=nwide)
  int    nwide;                    // number of wide band channels
  float  wfreq[MAXWIDE];           // freq
  float  wwidth[MAXWIDE];          // width
  // ifchain for ATCA
  int    chain[MAXWIN+MAXWIDE];    // if conversion chain


} WINDOW;

//  the maximum number of fields in mosaicing observations, careful, miriad
//  often has a smaller value, like 64 or 128.
#ifndef MAXFIELD
# define MAXFIELD  10000
#endif

#ifndef MAXMSG
# define MAXMSG 256
#endif

#ifndef MAXFSET
# define MAXFSET 100
#endif

// a helper class 

class Importmiriad
{
  // This is an implementation helper class used to store 'local' data
  // during the filling process.
public:
  // Create from a miriad dataset (a directory)
  Importmiriad(String& infile, Int debug=0, 
              Bool Qtsys=False,
              Bool Qarrays=False,
              Bool Qlinecal=False);

  // Standard destructor
  ~Importmiriad();
  
  // Check some of the contents of the data and header read
  void checkInput(Block<Int>& spw, Block<Int>& wide);

  // Debug output level
  Bool Debug(int level);

  // Set up the MeasurementSet, including StorageManagers and fixed columns.
  // If useTSM is True, the Tiled Storage Manager will be used to store
  // DATA, FLAG and WEIGHT_SPECTRUM
  void setupMeasurementSet(const String& MSFileName, Bool useTSM=True);

  // Fill the main table by reading in all the visibilities
  void fillMSMainTable();

  // Make an Antenna Table (can be called incrementally now)
  void fillAntennaTable();

  // Make a Syscal Table (can be called incrementally)
  void fillSyscalTable();

  // fill Spectralwindow table 
  void fillSpectralWindowTable(String vel);

  // fill Field table 
  void fillFieldTable();

  // fill Source table 
  void fillSourceTable();

  // fill the Feed table with minimal info needed for synthesis processing
  void fillFeedTable();

  // Fill the Observation and History (formerly ObsLog) tables
  void fillObsTables();

  // fix up the EPOCH MEASURE_REFERENCE keywords using the value found
  // in the (last) AN table - check if miriad really needs it
  void fixEpochReferences();

  void Tracking(int record);
  void check_window();
  Bool compareWindows(WINDOW& w1, WINDOW& w2);
  void Error(char *msg);
  void Warning(char *msg);
  void show();
  void close();

private:
  String                 infile_p;     // filename
  Int                    uv_handle_p;  // miriad handle 
  MeasurementSet         ms_p;         // the MS itself
  MSColumns             *msc_p;        // handy pointer to the columns in an MS  
  Int                    debug_p;      // debug level
  String                 array_p, 
                         project_p, 
                         object_p, 
                         telescope_p, 
                         observer_p, 
                         version_p,
                         timsys_p;
  Vector<Int>            nPixel_p, corrType_p, corrIndex_p;
  Matrix<Int>            corrProduct_p;
  Double                 epoch_p;
  MDirection::Types      epochRef_p;
  Int                    nArray_p;      // number of arrays (nAnt_p.nelements())
  Block<Int>             nAnt_p;        // number of antennas per array
  Block<Vector<Double> > receptorAngle_p;
  Vector<Double>         arrayXYZ_p;    // needs to be made with 3 elements
  Vector<Double>         ras_p, decs_p; // ra/dec for source list (source_p)
  Vector<String>         source_p,      // list of source names (?? object_p ??)
                         purpose_p;     // purpose of this source 
  LogIO                  os_p;          // logger


  // the following variables are for miriad, hence not Double/Int/Float
  // thus the code may have to be fixed on machines where these do not
  // agree ... may need special access code to get those into CASA
  // types on 64 bit machines??

  double preamble[5], first_time;
  int    ifield, nfield, npoint, nsource;     // both dra/ddec should become Vector's
  float  dra[MAXFIELD], ddec[MAXFIELD];       // offset in radians
  double ra[MAXFIELD], dec[MAXFIELD];
  int    field[MAXFIELD];                     // source index
  int    fcount[MAXFIELD];
  float  dra_p=0, ddec_p=0;
  int    pol_p[4];
  char   message[MAXMSG];


  // The following items more or less follow the uv variables in a dataset
  Int    nants_p, nants_offset_p, nchan_p, nwide_p, npol_p;
  Double antpos[3*MAXANT];
  Float  phasem1[MAXANT];
  Double ra_p, dec_p;       // current pointing center RA,DEC at EPOCH 
  Float  inttime_p, jyperk_p;
  Double freq_p;            // rest frequency of the primary line
  Int    mount_p;
  Double time_p;            // current MJD time
  Double timeFirst_p;       // First MJD time encountered

  // MIRIAD spectral window definition
  Int    freqSet_p,nFreqSet_p,ddid_p;
  WINDOW win[MAXFSET];  // allow for 16 different frequency setups
  
  Bool   Qtsys_p;    /* tsys weight's */
  Bool   Qarrays_p;  /* write separate arrays */
  Bool   Qlinecal_p; /* do linecal */
  Bool   keep[MAXWIN+MAXWIDE]; // keep this window for output to MS

  // Data buffers.... again in MIRIAD format
  
  float  data[2*MAXCHAN], wdata[2*MAXCHAN];     // 2*MAXCHAN since (Re,Im) pairs complex numbers
  int    flags[MAXCHAN], wflags[MAXCHAN];
  float  systemp[MAXANT*MAXWIDE];
  int    zero_tsys;

  // Counters
  int    nvis;
};
} // Namespace CASA - END

#endif