Line data Source code
1 : // -*- C++ -*-
2 : //# VisibilityResamplerBase.h: Definition of the VisibilityResamplerBase class
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 :
29 : #ifndef SYNTHESIS_VISIBILITYRESAMPLERBASE_H
30 : #define SYNTHESIS_VISIBILITYRESAMPLERBASE_H
31 :
32 : #include <synthesis/TransformMachines/CFStore.h>
33 : #include <synthesis/TransformMachines/CFStore2.h>
34 : #include <synthesis/TransformMachines/ConvolutionFunction.h>
35 : #include <synthesis/TransformMachines/Utils.h>
36 : #include <synthesis/TransformMachines/VBStore.h>
37 : #include <msvis/MSVis/VisBuffer.h>
38 : #include <casacore/casa/Arrays/Array.h>
39 : #include <casacore/casa/Arrays/Vector.h>
40 : #include <stdcasa/thread/AsynchronousTools.h>
41 :
42 : #include <casacore/casa/Logging/LogIO.h>
43 : #include <casacore/casa/Logging/LogSink.h>
44 : #include <casacore/casa/Logging/LogMessage.h>
45 : #include <casacore/casa/OS/Timer.h>
46 :
47 : namespace casa { //# NAMESPACE CASA - BEGIN
48 : class VisibilityResamplerBase
49 : {
50 : public:
51 0 : VisibilityResamplerBase():
52 0 : runTimeG_p(0.0), runTimeDG_p(0.0),runTimeG1_p(0.0), runTimeG2_p(0.0), runTimeG3_p(0.0), runTimeG4_p(0.0), runTimeG5_p(0.0), runTimeG6_p(0.0), runTimeG7_p(0.0),
53 0 : timer_p(),
54 0 : uvwScale_p(), offset_p(), chanMap_p(), polMap_p(), spwChanFreq_p(), spwChanConjFreq_p (), convFuncStore_p(), inc_p(),
55 0 : cfMap_p(), conjCFMap_p(), paTolerance_p(360.0)
56 :
57 0 : {};
58 : // VisibilityResamplerBase(const CFStore& cfs):
59 : // uvwScale_p(), offset_p(), chanMap_p(), polMap_p(), convFuncStore_p(), inc_p(),
60 : // cfMap_p(), conjCFMap_p()
61 : // {setConvFunc(cfs);};
62 :
63 : VisibilityResamplerBase(const VisibilityResamplerBase& other):
64 : uvwScale_p(), offset_p(), chanMap_p(), polMap_p(), spwChanFreq_p(), spwChanConjFreq_p (), convFuncStore_p(), inc_p(),
65 : cfMap_p(), conjCFMap_p(), paTolerance_p(360.0)
66 : {copy(other);}
67 :
68 0 : virtual ~VisibilityResamplerBase() {};
69 :
70 : VisibilityResamplerBase& operator=(const VisibilityResamplerBase& other);
71 :
72 : virtual VisibilityResamplerBase* clone() = 0;
73 :
74 : virtual void copy(const VisibilityResamplerBase& other);
75 : virtual void setParams(const casacore::Vector<casacore::Double>& uvwScale,
76 : const casacore::Vector<casacore::Double>& offset,
77 : const casacore::Vector<casacore::Double>& dphase) = 0;
78 :
79 : virtual void setMaps(const casacore::Vector<casacore::Int>& chanMap, const casacore::Vector<casacore::Int>& polMap) = 0;
80 : virtual void setCFMaps(const casacore::Vector<casacore::Int>& cfMap, const casacore::Vector<casacore::Int>& conjCFMap)=0;
81 : virtual void setFreqMaps(const casacore::Matrix<casacore::Double>& spwChanFreqs, const casacore::Matrix<casacore::Double>& spwnChanConjFreqs) = 0;
82 :
83 : virtual void setConvFunc(const CFStore& cfs) = 0;
84 : virtual void setPATolerance(const double& dPA) = 0;
85 : //
86 : //------------------------------------------------------------------------------
87 : //
88 : // Re-sample the griddedData on the VisBuffer (a.k.a gridding).
89 : //
90 : // In this class, these just call the private templated version.
91 : // The first variant grids onto a double precision grid while the
92 : // second one does it on a single precision grid.
93 : //
94 : virtual void DataToGrid(casacore::Array<casacore::DComplex>& griddedData, VBStore& vbs,
95 : casacore::Matrix<casacore::Double>& sumwt, const casacore::Bool& dopsf,
96 : casacore::Bool useConjFreqCF=false) = 0;
97 :
98 : virtual void DataToGrid(casacore::Array<casacore::Complex>& griddedData, VBStore& vbs,
99 : casacore::Matrix<casacore::Double>& sumwt, const casacore::Bool& dopsf,
100 : casacore::Bool useConjFreqCF=false) = 0;
101 : //
102 : //------------------------------------------------------------------------------
103 : //
104 : // Re-sample VisBuffer to a regular grid (griddedData) (a.k.a. de-gridding)
105 : //
106 : virtual void GridToData(VBStore& vbs,const casacore::Array<casacore::Complex>& griddedData) = 0;
107 : // virtual void GridToData(VBStore& vbs, casacore::Array<casacore::Complex>& griddedData);
108 :
109 : virtual void ComputeResiduals(VBStore& vbs) = 0;
110 :
111 : // Forward looking genealogical baggage -- required for the
112 : // MultiThreadedVisibilityResampler
113 : virtual void init(const casacore::Bool& doublePrecision) = 0;
114 : virtual void GatherGrids(casacore::Array<casacore::DComplex>& griddedData, casacore::Matrix<casacore::Double>& sumwt) = 0;
115 : virtual void GatherGrids(casacore::Array<casacore::Complex>& griddedData, casacore::Matrix<casacore::Double>& sumwt) = 0;
116 : virtual void initializePutBuffers(const casacore::Array<casacore::DComplex>& griddedData,
117 : const casacore::Matrix<casacore::Double>& sumwt) = 0;
118 : virtual void initializePutBuffers(const casacore::Array<casacore::Complex>& griddedData,
119 : const casacore::Matrix<casacore::Double>& sumwt) = 0;
120 : virtual void initializeDataBuffers(VBStore& vbs)=0;
121 : //
122 : // Aliases for more readable code at the FTMachine layer.
123 : //
124 0 : inline void finalizeToSky(casacore::Array<casacore::DComplex>& griddedData, casacore::Matrix<casacore::Double>& sumwt)
125 0 : {GatherGrids(griddedData, sumwt);};
126 0 : inline void finalizeToSky(casacore::Array<casacore::Complex>& griddedData, casacore::Matrix<casacore::Double>& sumwt)
127 0 : {GatherGrids(griddedData, sumwt);};
128 0 : inline void initializeToSky(const casacore::Array<casacore::DComplex>& griddedData,const casacore::Matrix<casacore::Double>& sumwt)
129 0 : {initializePutBuffers(griddedData, sumwt);};
130 0 : inline void initializeToSky(const casacore::Array<casacore::Complex>& griddedData,const casacore::Matrix<casacore::Double>& sumwt)
131 0 : {initializePutBuffers(griddedData, sumwt);};
132 : const casacore::Vector<casacore::Int> getCFMap() {return cfMap_p;};
133 : const casacore::Vector<casacore::Int> getConjCFMap() {return conjCFMap_p;};
134 :
135 :
136 : virtual void releaseBuffers() = 0;
137 : VBRow2CFMapType& getVBRow2CFMap() {return vbRow2CFMap_p;};
138 0 : VBRow2CFBMapType& getVBRow2CFBMap() {return vbRow2CFBMap_p;};
139 : virtual casacore::Int makeVBRow2CFMap(CFStore2& cfs,
140 : ConvolutionFunction& cf,
141 : const VisBuffer& vb, const casacore::Quantity& dPA,
142 : const casacore::Vector<casacore::Int>& dataChan2ImChanMap,
143 : const casacore::Vector<casacore::Int>& dataPol2ImPolMap,
144 : const casacore::Vector<casacore::Double>& pointingOffset);
145 :
146 : casacore::Double runTimeG_p, runTimeDG_p, runTimeG1_p, runTimeG2_p, runTimeG3_p, runTimeG4_p, runTimeG5_p, runTimeG6_p, runTimeG7_p;
147 : casacore::Timer timer_p;
148 : //
149 : //------------------------------------------------------------------------------
150 : //----------------------------Private parts-------------------------------------
151 : //------------------------------------------------------------------------------
152 : //
153 : // private:
154 : protected:
155 : casacore::Vector<casacore::Double> uvwScale_p, offset_p, dphase_p;
156 : casacore::Vector<casacore::Int> chanMap_p, polMap_p;
157 : casacore::Matrix<casacore::Double> spwChanFreq_p, spwChanConjFreq_p;
158 : CFStore convFuncStore_p;
159 : // casacore::Int inc0_p, inc1_p, inc2_p, inc3_p;
160 : casacore::Vector<casacore::Int> inc_p;
161 : casacore::Int* __restrict__ incPtr_p;
162 : casacore::Vector<casacore::Int> cfMap_p, conjCFMap_p;
163 : VBRow2CFMapType vbRow2CFMap_p;
164 : VBRow2CFBMapType vbRow2CFBMap_p;
165 : double paTolerance_p;
166 :
167 : void sgrid(casacore::Int& ndim,
168 : casacore::Double* __restrict__ pos,
169 : casacore::Int* __restrict__ loc,
170 : casacore::Int* __restrict__ off,
171 : casacore::Complex& phasor, const casacore::Int& irow,
172 : const casacore::Double* __restrict__ uvw,
173 : const casacore::Double& dphase, const casacore::Double& freq,
174 : const casacore::Double* __restrict__ scale,
175 : const casacore::Double* __restrict__ offset,
176 : const casacore::Float* __restrict__ sampling);
177 :
178 0 : inline casacore::Bool onGrid (const casacore::Int& nx, const casacore::Int& ny,
179 : const casacore::Vector<casacore::Int>& __restrict__ loc,
180 : const casacore::Vector<casacore::Int>& __restrict__ support) __restrict__
181 : {
182 0 : return (((loc(0)-support[0]) >= 0 ) && ((loc(0)+support[0]) < nx) &&
183 0 : ((loc(1)-support[1]) >= 0 ) && ((loc(1)+support[1]) < ny));
184 : };
185 : inline casacore::Bool onGrid (const casacore::Int& nx, const casacore::Int& ny,
186 : const casacore::Int& loc0, const casacore::Int& loc1,
187 : const casacore::Int& support) __restrict__
188 : {
189 : return (((loc0-support) >= 0 ) && ((loc0+support) < nx) &&
190 : ((loc1-support) >= 0 ) && ((loc1+support) < ny));
191 : };
192 :
193 : // casacore::Array assignment operator in CASACore requires lhs.nelements()
194 : // == 0 or lhs.nelements()=rhs.nelements()
195 : // template <class T>
196 : // inline void SETVEC(casacore::Vector<T>& lhs, const casacore::Vector<T>& rhs)
197 : // {lhs.resize(rhs.shape()); lhs = rhs;};
198 :
199 :
200 : //===============================================================================
201 : // CASACORE-LEVEL MATERIAL
202 : //===============================================================================
203 : // Internal methods to address a 4D array. These should ulimately
204 : // moved to a Array4D class in CASACore
205 : //
206 :
207 : // This is called less frequently. Currently once per VisBuffer
208 0 : inline void cacheAxisIncrements(const casacore::Int& n0, const casacore::Int& n1, const casacore::Int& n2, const casacore::Int& n3)
209 : {
210 : // inc0_p=1, inc1_p=inc0_p*n0, inc2_p=inc1_p*n1, inc3_p=inc2_p*n2;(void)n3;
211 0 : inc_p.resize(4);
212 0 : inc_p[0]=1; inc_p[1]=inc_p[0]*n0; inc_p[2]=inc_p[1]*n1; inc_p[3]=inc_p[2]*n2;(void)n3;
213 : casacore::Bool D;
214 0 : incPtr_p = inc_p.getStorage(D);
215 0 : }
216 0 : inline void cacheAxisIncrements(const casacore::Vector<casacore::Int>& n)
217 0 : {cacheAxisIncrements(n[0],n[1],n[2],n[3]);}
218 :
219 0 : inline void cacheAxisIncrements(const casacore::Vector<casacore::Int>& n, casacore::Vector<casacore::Int>& inc)
220 0 : {inc.resize(4);inc[0]=1; inc[1]=inc[0]*n[0]; inc[2]=inc[1]*n[1]; inc[3]=inc[2]*n[2];(void)n[3];}
221 :
222 : inline void cacheAxisIncrements(const casacore::Int n[4], casacore::Int inc[4])
223 : {inc[0]=1; inc[1]=inc[0]*n[0]; inc[2]=inc[1]*n[1]; inc[3]=inc[2]*n[2];(void)n[3];}
224 :
225 : // Version that use internally cached inc_p
226 : // template <class T>
227 0 : inline void addTo4DArray(casacore::DComplex* __restrict__& store, casacore::Int* __restrict__& iPos,
228 : casacore::Complex& nvalue, casacore::Double& wt) __restrict__
229 0 : {addTo4DArray(store, iPos, incPtr_p, nvalue, wt);}
230 :
231 0 : inline void addTo4DArray(casacore::Complex* __restrict__& store, casacore::Int* __restrict__& iPos,
232 : casacore::Complex& nvalue, casacore::Double& wt) __restrict__
233 0 : {addTo4DArray(store, iPos, incPtr_p, nvalue, wt);}
234 :
235 :
236 : // Version where inc_p is supplied from outside
237 0 : inline void addTo4DArray(casacore::DComplex* __restrict__& store, casacore::Int* __restrict__& iPos,
238 : casacore::Int* __restrict__ inc, casacore::Complex& nvalue, casacore::Double& wt) __restrict__
239 0 : {store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]] += (nvalue*casacore::Complex(wt));}
240 :
241 0 : inline void addTo4DArray(casacore::Complex* __restrict__& store, casacore::Int* __restrict__& iPos,
242 : casacore::Int* __restrict__ inc, casacore::Complex& nvalue, casacore::Double& wt) __restrict__
243 0 : {store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]] += (nvalue*casacore::Complex(wt));}
244 :
245 :
246 0 : inline casacore::Complex getFrom4DArray(const casacore::Complex* __restrict__& store,
247 : const casacore::Int* __restrict__& iPos,
248 : const casacore::Vector<casacore::Int>& inc)
249 : // __restrict__
250 0 : {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
251 :
252 : inline casacore::Complex getFrom4DArray(const casacore::Complex* __restrict__& store,
253 : const casacore::Vector<casacore::Int> iPos, const casacore::Vector<casacore::Int>& inc)
254 : // __restrict__
255 : {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
256 :
257 : inline casacore::DComplex getFrom4DArray(const casacore::DComplex* __restrict__& store,
258 : const casacore::Int* __restrict__& iPos,
259 : const casacore::Vector<casacore::Int>& inc)
260 : // __restrict__
261 : {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
262 :
263 : inline casacore::DComplex getFrom4DArray(const casacore::DComplex* __restrict__& store,
264 : const casacore::Vector<casacore::Int> iPos, const casacore::Vector<casacore::Int>& inc)
265 : // __restrict__
266 : {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
267 :
268 :
269 : // The following two methods are called in the innermost loop.
270 0 : inline casacore::Complex getFrom4DArray(const casacore::Complex* __restrict__& store, const casacore::Int* __restrict__& iPos)
271 : // __restrict__
272 0 : {return getFrom4DArray(store, iPos, inc_p);}
273 :
274 : inline casacore::DComplex getFrom4DArray(const casacore::DComplex* __restrict__& store, const casacore::Int* __restrict__& iPos)
275 : // __restrict__
276 : {return getFrom4DArray(store, iPos, inc_p);}
277 :
278 : };
279 : }; //# NAMESPACE CASA - END
280 :
281 : #endif //
|
