Line data Source code
1 : // -*- C++ -*-
2 : //# Utils.h: Definition of global functions in Utils.cc
3 : //# Copyright (C) 1997,1998,1999,2000,2001,2002,2003
4 : //# Associated Universities, Inc. Washington DC, USA.
5 : //#
6 : //# This library is free software; you can redistribute it and/or modify it
7 : //# under the terms of the GNU Library General Public License as published by
8 : //# the Free Software Foundation; either version 2 of the License, or (at your
9 : //# option) any later version.
10 : //#
11 : //# This library is distributed in the hope that it will be useful, but WITHOUT
12 : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 : //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
14 : //# License for more details.
15 : //#
16 : //# You should have received a copy of the GNU Library General Public License
17 : //# along with this library; if not, write to the Free Software Foundation,
18 : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
19 : //#
20 : //# Correspondence concerning AIPS++ should be addressed as follows:
21 : //# Internet email: casa-feedback@nrao.edu.
22 : //# Postal address: AIPS++ Project Office
23 : //# National Radio Astronomy Observatory
24 : //# 520 Edgemont Road
25 : //# Charlottesville, VA 22903-2475 USA
26 : //#
27 : //# $Id$
28 : #ifndef SYNTHESIS_TRANSFORM2_UTILS_H
29 : #define SYNTHESIS_TRANSFORM2_UTILS_H
30 :
31 :
32 : #include <casacore/casa/aips.h>
33 : #include <casacore/casa/Exceptions/Error.h>
34 : #include <msvis/MSVis/VisBuffer2.h>
35 :
36 : #include <casacore/images/Images/ImageOpener.h>
37 : #include <casacore/casa/Quanta/Quantum.h>
38 : #include <casacore/images/Images/ImageInterface.h>
39 : //#include <ms/MeasurementSets/MeasurementSet.h>
40 : #include <msvis/MSVis/VisibilityIterator2.h>
41 : #include <casacore/ms/MeasurementSets/MSColumns.h>
42 : #include <synthesis/TransformMachines/CFCell.h>
43 : #include <casacore/images/Images/TempImage.h>
44 : #include <casacore/casa/Arrays/Array.h>
45 : #include <casacore/casa/Logging/LogIO.h>
46 : #include <iostream>
47 : #include <array>
48 :
49 : namespace casa
50 : {
51 : using namespace vi;
52 : namespace refim {
53 : casacore::Int getPhaseCenter(casacore::MeasurementSet& ms, casacore::MDirection& dir0, casacore::Int whichField=-1);
54 : casacore::Bool findMaxAbsLattice(const casacore::ImageInterface<casacore::Float>& lattice,
55 : casacore::Float& maxAbs,casacore::IPosition& posMaxAbs);
56 : casacore::Bool findMaxAbsLattice(const casacore::ImageInterface<casacore::Float>& masklat,
57 : const casacore::Lattice<casacore::Float>& lattice,
58 : casacore::Float& maxAbs,casacore::IPosition& posMaxAbs,
59 : casacore::Bool flip=false);
60 : casacore::Double getCurrentTimeStamp(const VisBuffer2& vb);
61 : void makeStokesAxis(casacore::Int npol_p, casacore::Vector<casacore::String>& polType, casacore::Vector<casacore::Int>& whichStokes);
62 : casacore::Double getPA(const vi::VisBuffer2& vb);
63 : void storeImg(casacore::String fileName,casacore::ImageInterface<casacore::Complex>& theImg, casacore::Bool writeReIm=false);
64 : void storeImg(casacore::String fileName,casacore::ImageInterface<casacore::Float>& theImg);
65 : void storeArrayAsImage(casacore::String fileName, const casacore::CoordinateSystem& coords, const casacore::Array<casacore::Complex>& cf);
66 : void storeArrayAsImage(casacore::String fileName, const casacore::CoordinateSystem& coords, const casacore::Array<casacore::DComplex>& cf);
67 : void storeArrayAsImage(casacore::String fileName, const casacore::CoordinateSystem& coords, const casacore::Array<casacore::Float>& cf);
68 : void storeArrayAsImage(casacore::String fileName, const casacore::CoordinateSystem& coords, const casacore::Array<casacore::Double>& cf);
69 :
70 : casacore::Bool isVBNaN(const VisBuffer2& vb, casacore::String& mesg);
71 : namespace SynthesisUtils
72 : {
73 : //using namespace vi;
74 : void rotateComplexArray(casacore::LogIO& logIO, casacore::Array<casacore::Complex>& inArray,
75 : casacore::CoordinateSystem& inCS,
76 : casacore::Array<casacore::Complex>& outArray,
77 : casacore::Double dAngleRad,
78 : casacore::String interpMathod=casacore::String("CUBIC"),
79 : casacore::Bool modifyInCS=true);
80 : void findLatticeMax(const casacore::Array<casacore::Complex>& lattice,
81 : casacore::Vector<casacore::Float>& maxAbs,
82 : casacore::Vector<casacore::IPosition>& posMaxAbs) ;
83 : void findLatticeMax(const casacore::ImageInterface<casacore::Complex>& lattice,
84 : casacore::Vector<casacore::Float>& maxAbs,
85 : casacore::Vector<casacore::IPosition>& posMaxAbs) ;
86 : void findLatticeMax(const casacore::ImageInterface<casacore::Float>& lattice,
87 : casacore::Vector<casacore::Float>& maxAbs,
88 : casacore::Vector<casacore::IPosition>& posMaxAbs) ;
89 198448534 : inline casacore::Int nint(const casacore::Double& v) {return (casacore::Int)std::floor(v+0.5);}
90 3767304 : inline casacore::Int nint(const casacore::Float& v) {return (casacore::Int)std::floor(v+0.5);}
91 : inline casacore::Bool near(const casacore::Double& d1, const casacore::Double& d2,
92 : const casacore::Double EPS=1E-6)
93 : {
94 : casacore::Bool b1=(fabs(d1-d2) < EPS)?true:false;
95 : return b1;
96 : }
97 : template <class T>
98 50463 : inline void SETVEC(casacore::Vector<T>& lhs, const casacore::Vector<T>& rhs)
99 50463 : {lhs.resize(rhs.shape()); lhs = rhs;};
100 : template <class T>
101 14418 : inline void SETVEC(casacore::Array<T>& lhs, const casacore::Array<T>& rhs)
102 14418 : {lhs.resize(rhs.shape()); lhs = rhs;};
103 :
104 : template <class T>
105 : T getenv(const char *name, const T defaultVal);
106 : casacore::Float libreSpheroidal(casacore::Float nu);
107 : casacore::Double getRefFreq(const VisBuffer2& vb);
108 : void makeFTCoordSys(const casacore::CoordinateSystem& coords,
109 : const casacore::Int& convSize,
110 : const casacore::Vector<casacore::Double>& ftRef,
111 : casacore::CoordinateSystem& ftCoords);
112 :
113 : void expandFreqSelection(const casacore::Matrix<casacore::Double>& freqSelection,
114 : casacore::Matrix<casacore::Double>& expandedFreqList,
115 : casacore::Matrix<casacore::Double>& expandedConjFreqList);
116 :
117 : template <class T>
118 : void libreConvolver(casacore::Array<T>& c1, const casacore::Array<T>& c2);
119 591 : inline static casacore::Double conjFreq(const casacore::Double& freq, const casacore::Double& refFreq)
120 591 : {return sqrt(2*refFreq*refFreq - freq*freq);};
121 :
122 : casacore::Double nearestValue(const casacore::Vector<casacore::Double>& list, const casacore::Double& val, casacore::Int& index);
123 :
124 : template <class T>
125 : T stdNearestValue(const std::vector<T>& list, const T& val, casacore::Int& index);
126 :
127 : casacore::CoordinateSystem makeUVCoords(const casacore::CoordinateSystem& imageCoordSys,
128 : const casacore::IPosition& shape);
129 :
130 : casacore::Vector<casacore::Int> mapSpwIDToDDID(const VisBuffer2& vb, const casacore::Int& spwID);
131 : casacore::Vector<casacore::Int> mapSpwIDToPolID(const VisBuffer2& vb, const casacore::Int& spwID);
132 : void calcIntersection(const casacore::Int blc1[2], const casacore::Int trc1[2], const casacore::Float blc2[2], const casacore::Float trc2[2],
133 : casacore::Float blc[2], casacore::Float trc[2]);
134 : casacore::Bool checkIntersection(const casacore::Int blc1[2], const casacore::Int trc1[2], const casacore::Float blc2[2], const casacore::Float trc2[2]);
135 :
136 : casacore::String mjdToString(casacore::Time& mjd);
137 :
138 : template<class Iterator>
139 : Iterator Unique(Iterator first, Iterator last);
140 :
141 : void showCS(const casacore::CoordinateSystem& cs, std::ostream& os, const casacore::String& msg=casacore::String());
142 : const casacore::Array<casacore::Complex> getCFPixels(const casacore::String& Dir, const casacore::String& fileName);
143 : void putCFPixels(const casacore::String& Dir, const casacore::String& fileName,
144 : const casacore::Array<casacore::Complex>& srcpix);
145 : const casacore::IPosition getCFShape(const casacore::String& Dir, const casacore::String& fileName);
146 :
147 : void rotate2(const double& actualPA, CFCell& baseCFC, CFCell& cfc, const double& rotAngleIncr);
148 :
149 : casacore::TableRecord getCFParams(const casacore::String& dirName,const casacore::String& fileName,
150 : casacore::IPosition& cfShape,
151 : casacore::Array<casacore::Complex>& pixelBuffer,
152 : casacore::CoordinateSystem& coordSys,
153 : casacore::Double& sampling,
154 : casacore::Double& paVal,
155 : casacore::Int& xSupport, casacore::Int& ySupport,
156 : casacore::Double& fVal, casacore::Double& wVal, casacore::Int& mVal,
157 : casacore::Double& conjFreq, casacore::Int& conjPoln,
158 : casacore::Bool loadPixels,
159 : casacore::Bool loadMiscInfo=true);
160 :
161 : casacore::Vector<casacore::String> parseBandName(const casacore::String& fullName);
162 :
163 : casacore::CoordinateSystem makeModelGridFromImage(const std::string& modelImageName,
164 : casacore::TempImage<casacore::DComplex>& modelImageGrid);
165 :
166 : void makeAWLists(const casacore::Vector<double>& wVals,
167 : const casacore::Vector<double>& fVals,
168 : const bool& wbAWP, const uint& nw,
169 : const double& imRefFreq, const double& spwRefFreq,
170 : casacore::Vector<int>& wNdxList, casacore::Vector<int>& spwNdxList,
171 : const int vbSPW);
172 :
173 :
174 :
175 :
176 :
177 : }; ///end of namespace SynthesisUtils
178 :
179 : class MathUtils
180 : {
181 : public:
182 : MathUtils();
183 : casacore::Float interpLanczos( const casacore::Double& x , const casacore::Double& y, const casacore::Double& nx, const casacore::Double& ny, const casacore::Float* data, const casacore::Float a=3);
184 : casacore::Float sinc(const casacore::Float x) ;
185 : casacore::Array<casacore::Complex> resample(const casacore::Array<casacore::Complex>& inarray, const casacore::Double factorX, const casacore::Double factorY);
186 : static casacore::Array<casacore::Complex> resampleViaFFT(const casacore::Array<casacore::Complex>& inarray, const casacore::Double factorX, const casacore::Double factorY);
187 : static casacore::Array<casacore::Complex> resampleViaFFT(const casacore::Array<casacore::Complex>& inarray, const casacore::Int nX, const casacore::Int nY);
188 : //get the middle nx, ny of inArr ...obviously inArr fist 2 dim must be bigger or
189 : //equal to nx ny
190 : static casacore::Array<casacore::Complex> getMiddle(const casacore::Array<casacore::Complex>& inArr, const int nx, const int ny);
191 : //put inArr in the middle of ouArr...first 2 dim of inArr have to be smaller or equal to outArr and all other axes have to be of the same length
192 : static void putMiddle(casacore::Array<casacore::Complex>& outArr, const casacore::Array<casacore::Complex>& inArr);
193 : private:
194 : std::array<casacore::Float,8000> sincCache_p;
195 : casacore::Float *sincCachePtr_p;
196 : void initSincCache();
197 :
198 :
199 :
200 : };
201 :
202 : void getHADec(casacore::MeasurementSet& ms, const VisBuffer2& vb, casacore::Double &HA, casacore::Double& RA, casacore::Double& Dec);
203 :
204 : /////////////////////////////////////////////////////////////////////////////
205 : //
206 : // An interface class to detect changes in the VisBuffer
207 : // Exact meaning of the "change" is defined in the derived classes
208 : //
209 : struct IChangeDetector {
210 : // return true if a change occurs in the given row since the last call of update
211 : virtual casacore::Bool changed(const VisBuffer2 &vb, casacore::Int row) const = 0;
212 : // start looking for a change from the given row of the VisBuffer
213 : virtual void update(const VisBuffer2 &vb, casacore::Int row) = 0;
214 :
215 : // reset to the state which exists just after construction
216 : virtual void reset() = 0;
217 :
218 : // some derived methods, which use the abstract virtual function changed(vb,row)
219 :
220 : // return true if a change occurs somewhere in the buffer
221 : casacore::Bool changed(const VisBuffer2 &vb) const;
222 : // return true if a change occurs somewhere in the buffer starting from row1
223 : // up to row2 (row2=-1 means up to the end of the buffer). The row number,
224 : // where the change occurs is returned in the row2 parameter
225 : casacore::Bool changedBuffer(const VisBuffer2 &vb, casacore::Int row1, casacore::Int &row2) const;
226 : protected:
227 : // a virtual destructor to make the compiler happy
228 : virtual ~IChangeDetector();
229 : };
230 : //
231 : //////////////////////////////////////////////////////////////////////////////
232 :
233 : //////////////////////////////////////////////////////////////////////////////
234 : //
235 : // ParAngleChangeDetector - a class to detect a change in the parallactic
236 : // angle.
237 : //
238 : class ParAngleChangeDetector : public IChangeDetector {
239 : casacore::Double pa_tolerance_p; // a parallactic angle tolerance. If exeeded,
240 : // the angle is considered to be changed.
241 : casacore::Double last_pa_p; // last value of the parallactic angle
242 : public:
243 : // The default constructor
244 352 : ParAngleChangeDetector():pa_tolerance_p(0.0) {};
245 : // set up the tolerance, which determines how much the position angle should
246 : // change to report the change by this class
247 : ParAngleChangeDetector(const casacore::Quantity &pa_tolerance);
248 :
249 : virtual void setTolerance(const casacore::Quantity &pa_tolerance);
250 : // reset to the state which exists just after construction
251 : virtual void reset();
252 :
253 : // return parallactic angle tolerance
254 : casacore::Quantity getParAngleTolerance() const;
255 :
256 : // implementation of the base class' virtual functions
257 :
258 : // return true if a change occurs in the given row since the last call of update
259 : virtual casacore::Bool changed(const VisBuffer2 &vb, casacore::Int row) const;
260 : // start looking for a change from the given row of the VisBuffer
261 : virtual void update(const VisBuffer2 &vb, casacore::Int row);
262 : };
263 :
264 : //
265 : /////////////////////////////////////////////////////////////////////////////
266 :
267 : };
268 : };
269 : #endif
|