casa / casa5 / 2267d148f27

Commits

Ville Suoranta authored 2267d148f27 Merge 21 Jun 2017

Merge pull request #118 in CASA/casa from bugfix/CAS-10280 to master

* commit '80489706b417542a1533129516e13408111575da':
  Sped up the define image by passing the channel info to SynthesisUtilsMethod
  changed the code for some frequency range calculation not to use the iterator but use MSUtil functionalities

code/msvis/MSVis/MSUtil.cc

Modified

-172
+				  outnchan[nout]=locnchan[j];
-173
-174
-175
-176
+			  }
-177
+		  }
-178
-179
+	  }
-180
-181
-+
+  }
-+
+  void MSUtil::getChannelRangeFromFreqRange(Int& start,
-+
+				  Int& nchan,
-+
+				  const MeasurementSet& ms,
-+
+				  const Int spw,
-+
+				  const Double freqStart,
-+
+				  const Double freqEnd,
-+
+				  const Double freqStep,
-+
+			    const MFrequency::Types freqframe)
-+
+  {
-+
+    start=-1;
-+
+    nchan=-1;
-+
+    ROMSFieldColumns fieldCol(ms.field());
-+
+    ROMSDataDescColumns ddCol(ms.dataDescription());
-+
+	Vector<Int> dataDescSel=MSDataDescIndex(ms.dataDescription()).matchSpwId(spw);
-+
+	//cerr << "dataDescSel " << dataDescSel << endl;
-+
+	if(dataDescSel.nelements()==0)
-+
+		return;
-+
+	Vector<Double> t;
-+
+    ROScalarColumn<Double> (ms,MS::columnName(MS::TIME)).getColumn(t);
-+
+	Vector<Int> ddId;
-+
+    ROScalarColumn<Int> (ms,MS::columnName(MS::DATA_DESC_ID)).getColumn(ddId);
-+
+    ROMSSpWindowColumns spwCol(ms.spectralWindow());
-+
+    Vector<Double> ddIdD(ddId.shape());
-+
+    convertArray(ddIdD, ddId);
-+
+    ddIdD+= 1.0; //no zero id
-+
+     //we have to do this as one can have multiple dd for the same time.
-+
+    t*=ddIdD;
-+
+    //t(fldId != fieldId)=-1.0;
-+
+    Vector<Double> elt;
-+
+    Vector<Int> elindx;
-+
+    //rejecting the large blocks of same time for all baselines
-+
+    //this speeds up by a lot GenSort::sort
-+
+    rejectConsecutive(t, elt, elindx);
-+
-+
+    ROScalarMeasColumn<MEpoch> timeCol(ms, MS::columnName(MS::TIME));
-+
+    Vector<uInt>  uniqIndx;
-+
+    uInt nTimes=GenSortIndirect<Double>::sort (uniqIndx, elt, Sort::Ascending, Sort::QuickSort|Sort::NoDuplicates);
-+
-+
+    t.resize(0);
-+
-+
+    //ROScalarColumn<Int> ddId(ms,MS::columnName(MS::DATA_DESC_ID));
-+
+    ROScalarColumn<Int> fldId(ms,MS::columnName(MS::FIELD_ID));
-+
+    //now need to do the conversion to data frame from requested frame
-+
+    //Get the epoch mesasures of the first row
-+
+    MEpoch ep;
-+
+    timeCol.get(0, ep);
-+
+    String observatory;
-+
+    MPosition obsPos;
-+
+    /////observatory position
-+
+    ROMSColumns msc(ms);
-+
+    if (ms.observation().nrow() > 0) {
-+
+      observatory = msc.observation().telescopeName()(msc.observationId()(0));
-+
+    }
-+
+    if (observatory.length() == 0 ||
-+
+	!MeasTable::Observatory(obsPos,observatory)) {
-+
+      // unknown observatory, use first antenna
-+
+      obsPos=msc.antenna().positionMeas()(0);
-+
+    }
-+
+    //////
-+
+    Int oldDD=dataDescSel(0);
-+
+    Int newDD=oldDD;
-+
+    //For now we will assume that the field is not moving very far from polynome 0
-+
+    MDirection dir =fieldCol.phaseDirMeas(0);
-+
+    MFrequency::Types obsMFreqType= (MFrequency::Types) (spwCol.measFreqRef()(ddCol.spectralWindowId()(dataDescSel(0))));
-+
+    MeasFrame frame(ep, obsPos, dir);
-+
+    MFrequency::Convert toObs(freqframe,
-+
+                              MFrequency::Ref(obsMFreqType, frame));
-+
+    Double freqEndMax=freqEnd+0.5*fabs(freqStep);
-+
+    Double freqStartMin=freqStart-0.5*fabs(freqStep);
-+
+    if(freqframe != obsMFreqType){
-+
+      freqEndMax=0.0;
-+
+      freqStartMin=C::dbl_max;
-+
+    }
-+
+    for (uInt j=0; j< nTimes; ++j) {
-+
+        if(anyEQ(dataDescSel, ddId(elindx[uniqIndx[j]]))){
-+
+            timeCol.get(elindx[uniqIndx[j]], ep);
-+
+            newDD=ddId(elindx[uniqIndx[j]]);
-+
+            if(oldDD != newDD) {
-+
-+
+                oldDD=newDD;
-+
+                if(spwCol.measFreqRef()(ddCol.spectralWindowId()(newDD)) != obsMFreqType) {
-+
+                    obsMFreqType= (MFrequency::Types) (spwCol.measFreqRef()(ddCol.spectralWindowId()(newDD)));
-+
+                    toObs.setOut(MFrequency::Ref(obsMFreqType, frame));
-+
+                }
-+
+            }
-+
+            if(obsMFreqType != freqframe) {
-+
+                dir=fieldCol.phaseDirMeas(fldId(elindx[uniqIndx[j]]));
-+
+                frame.resetEpoch(ep);
-+
+                frame.resetDirection(dir);
-+
+                Double freqTmp=toObs(Quantity(freqStart, "Hz")).get("Hz").getValue();
-+
+                freqStartMin=(freqStartMin > freqTmp) ? freqTmp : freqStartMin;
-+
+                freqTmp=toObs(Quantity(freqEnd, "Hz")).get("Hz").getValue();
-+
+                freqEndMax=(freqEndMax < freqTmp) ? freqTmp : freqEndMax;
-+
+            }
-+
+        }
-+
+    }
-+
-+
-+
+    MSSpwIndex spwIn(ms.spectralWindow());
-+
+	Vector<Int> spws;
-+
+	Vector<Int> starts;
-+
+	Vector<Int> nchans;
-+
+    if(!spwIn.matchFrequencyRange(freqStartMin-0.5*freqStep, freqEndMax+0.5*freqStep, spws, starts, nchans)){
-+
+			return;
-+
+	}
-+
+	for (uInt k=0; k < spws.nelements(); ++k){
-+
+			if(spws[k]==spw){
-+
+					start=starts[k];
-+
+					nchan=nchans[k];
-+
+			}
-+
-+
+	}
-295
+  }
-296
+  void MSUtil::getFreqRangeInSpw( Double& freqStart,
-297
+				  Double& freqEnd, const Vector<Int>& spw, const Vector<Int>& start,
-298
+				  const Vector<Int>& nchan,
-299
+				  const MeasurementSet& ms,
-300
+				  const MFrequency::Types freqframe,
--
+				  const Int fieldId){
-+
+				  const Int fieldId, const Bool fromEdge){
-302
-303
-304
+    freqStart=C::dbl_max;
-305
+    freqEnd=0.0;
-306
+    Vector<Double> t;
-307
+    ROScalarColumn<Double> (ms,MS::columnName(MS::TIME)).getColumn(t);
-308
+    Vector<Int> ddId;
-309
+    Vector<Int> fldId;
-310
+    ROScalarColumn<Int> (ms,MS::columnName(MS::DATA_DESC_ID)).getColumn(ddId);
-311
+    ROScalarColumn<Int> (ms,MS::columnName(MS::FIELD_ID)).getColumn(fldId);
-342
+    if (observatory.length() == 0 ||
-343
+	!MeasTable::Observatory(obsPos,observatory)) {
-344
+      // unknown observatory, use first antenna
-345
+      obsPos=msc.antenna().positionMeas()(0);
-346
+    }
-347
+    //////
-348
+    MeasFrame frame(ep, obsPos, dir);
-349
-350
-351
+    for (uInt ispw =0 ; ispw < spw.nelements() ; ++ispw){
-+
+		if(nchan[ispw]>0 && start[ispw] >-1){
-353
+      Double freqStartObs=C::dbl_max;
-354
+      Double freqEndObs=0.0;
-355
+      Vector<Double> chanfreq=spwCol.chanFreq()(spw[ispw]);
-356
+      Vector<Double> chanwid=spwCol.chanWidth()(spw[ispw]);
-357
+      Vector<Int> ddOfSpw=mddin.matchSpwId(spw[ispw]);
-358
+      for (Int ichan=start[ispw]; ichan<start[ispw]+nchan[ispw]; ++ichan){
--
+	if(freqStartObs > (chanfreq[ichan]-fabs(chanwid[ichan])/2.0)) freqStartObs=chanfreq[ichan]-fabs(chanwid[ichan])/2.0;
--
+	if(freqEndObs < (chanfreq[ichan]+fabs(chanwid[ichan])/2.0)) freqEndObs=chanfreq[ichan]+fabs(chanwid[ichan])/2.0;
-+
+		  if(fromEdge){
-+
+			if(freqStartObs > (chanfreq[ichan]-fabs(chanwid[ichan])/2.0)) freqStartObs=chanfreq[ichan]-fabs(chanwid[ichan])/2.0;
-+
+			if(freqEndObs < (chanfreq[ichan]+fabs(chanwid[ichan])/2.0)) freqEndObs=chanfreq[ichan]+fabs(chanwid[ichan])/2.0;
-+
+		  }
-+
+		  else{
-+
+			if(freqStartObs > (chanfreq[ichan])) freqStartObs=chanfreq[ichan];
-+
+			if(freqEndObs < (chanfreq[ichan])) freqEndObs=chanfreq[ichan];
-+
+		  }
-367
+      }
-368
+      obsMFreqType= (MFrequency::Types) (spwCol.measFreqRef()(spw[ispw]));
-369
+      if((obsMFreqType==MFrequency::REST) || (obsMFreqType==freqframe && obsMFreqType != MFrequency::TOPO)){
-370
+	if(freqStart > freqStartObs)  freqStart=freqStartObs;
-371
+	if(freqEnd < freqStartObs)  freqEnd=freqStartObs;
-372
+	if(freqStart > freqEndObs)  freqStart=freqEndObs;
-373
+	if(freqEnd < freqEndObs)  freqEnd=freqEndObs;
-374
+      }
-375
+      else{
-376
+	MFrequency::Convert toframe(obsMFreqType,
-382
+	    Double freqTmp=toframe(Quantity(freqStartObs, "Hz")).get("Hz").getValue();
-383
+	    if(freqStart > freqTmp)  freqStart=freqTmp;
-384
+	    if(freqEnd < freqTmp)  freqEnd=freqTmp;
-385
+	    freqTmp=toframe(Quantity(freqEndObs, "Hz")).get("Hz").getValue();
-386
+	    if(freqStart > freqTmp)  freqStart=freqTmp;
-387
+	    if(freqEnd < freqTmp)  freqEnd=freqTmp;
-388
+	  }
-389
+	}
-390
+      }
-391
+    }
-+
+	}
-393
+  }
-394
+  void MSUtil::rejectConsecutive(const Vector<Double>& t, Vector<Double>& retval, Vector<Int>& indx){
-395
+    uInt n=t.nelements();
-396
+    if(n >0){
-397
+      retval.resize(n);
-398
+      indx.resize(n);
-399
+      retval[0]=t[0];
-400
+      indx[0]=0;
-401
+    }
-402
+    else

casa5

Commits

Jira Issues

Add shortcut