00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029 #ifndef SYNTHESIS_TRANSFORM2_VISIBILITYRESAMPLERBASE_H
00030 #define SYNTHESIS_TRANSFORM2_VISIBILITYRESAMPLERBASE_H
00031
00032 #include <synthesis/TransformMachines2/CFStore.h>
00033 #include <synthesis/TransformMachines2/CFStore2.h>
00034 #include <synthesis/TransformMachines2/ConvolutionFunction.h>
00035 #include <synthesis/TransformMachines2/Utils.h>
00036 #include <synthesis/TransformMachines2/VBStore.h>
00037 #include <msvis/MSVis/VisBuffer2.h>
00038 #include <casa/Arrays/Array.h>
00039 #include <casa/Arrays/Vector.h>
00040 #include <msvis/MSVis/AsynchronousTools.h>
00041
00042 #include <casa/Logging/LogIO.h>
00043 #include <casa/Logging/LogSink.h>
00044 #include <casa/Logging/LogMessage.h>
00045 #include <casa/OS/Timer.h>
00046
00047 namespace casa {
00048 using namespace vi;
00049 namespace refim{
00050 class VisibilityResamplerBase
00051 {
00052 public:
00053 VisibilityResamplerBase():
00054 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),
00055 timer_p(),
00056 uvwScale_p(), offset_p(), chanMap_p(), polMap_p(), spwChanFreq_p(), spwChanConjFreq_p (), convFuncStore_p(), inc_p(),
00057 cfMap_p(), conjCFMap_p()
00058
00059 {};
00060
00061
00062
00063
00064
00065 VisibilityResamplerBase(const VisibilityResamplerBase& other):
00066 uvwScale_p(), offset_p(), chanMap_p(), polMap_p(), spwChanFreq_p(), spwChanConjFreq_p (), convFuncStore_p(), inc_p(),
00067 cfMap_p(), conjCFMap_p()
00068 {copy(other);}
00069
00070 virtual ~VisibilityResamplerBase() {};
00071
00072 VisibilityResamplerBase& operator=(const VisibilityResamplerBase& other);
00073
00074 virtual VisibilityResamplerBase* clone() = 0;
00075
00076 virtual void copy(const VisibilityResamplerBase& other);
00077 virtual void setParams(const Vector<Double>& uvwScale,
00078 const Vector<Double>& offset,
00079 const Vector<Double>& dphase) = 0;
00080
00081 virtual void setMaps(const Vector<Int>& chanMap, const Vector<Int>& polMap) = 0;
00082 virtual void setCFMaps(const Vector<Int>& cfMap, const Vector<Int>& conjCFMap)=0;
00083 virtual void setFreqMaps(const Matrix<Double>& spwChanFreqs, const Matrix<Double>& spwnChanConjFreqs) = 0;
00084
00085 virtual void setConvFunc(const CFStore& cfs) = 0;
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095 virtual void DataToGrid(Array<DComplex>& griddedData, VBStore& vbs,
00096 Matrix<Double>& sumwt, const Bool& dopsf,
00097 Bool useConjFreqCF=False) = 0;
00098
00099 virtual void DataToGrid(Array<Complex>& griddedData, VBStore& vbs,
00100 Matrix<Double>& sumwt, const Bool& dopsf,
00101 Bool useConjFreqCF=False) = 0;
00102
00103
00104
00105
00106
00107 virtual void GridToData(VBStore& vbs,const Array<Complex>& griddedData) = 0;
00108
00109
00110 virtual void ComputeResiduals(VBStore& vbs) = 0;
00111
00112
00113
00114 virtual void init(const Bool& doublePrecision) = 0;
00115 virtual void GatherGrids(Array<DComplex>& griddedData, Matrix<Double>& sumwt) = 0;
00116 virtual void GatherGrids(Array<Complex>& griddedData, Matrix<Double>& sumwt) = 0;
00117 virtual void initializePutBuffers(const Array<DComplex>& griddedData,
00118 const Matrix<Double>& sumwt) = 0;
00119 virtual void initializePutBuffers(const Array<Complex>& griddedData,
00120 const Matrix<Double>& sumwt) = 0;
00121 virtual void initializeDataBuffers(VBStore& vbs)=0;
00122
00123
00124
00125 inline void finalizeToSky(Array<DComplex>& griddedData, Matrix<Double>& sumwt)
00126 {GatherGrids(griddedData, sumwt);};
00127 inline void finalizeToSky(Array<Complex>& griddedData, Matrix<Double>& sumwt)
00128 {GatherGrids(griddedData, sumwt);};
00129 inline void initializeToSky(const Array<DComplex>& griddedData,const Matrix<Double>& sumwt)
00130 {initializePutBuffers(griddedData, sumwt);};
00131 inline void initializeToSky(const Array<Complex>& griddedData,const Matrix<Double>& sumwt)
00132 {initializePutBuffers(griddedData, sumwt);};
00133 const Vector<Int> getCFMap() {return cfMap_p;};
00134 const Vector<Int> getConjCFMap() {return conjCFMap_p;};
00135
00136
00137 virtual void releaseBuffers() = 0;
00138 VBRow2CFMapType& getVBRow2CFMap() {return vbRow2CFMap_p;};
00139 VBRow2CFBMapType& getVBRow2CFBMap() {return vbRow2CFBMap_p;};
00140 virtual Int makeVBRow2CFMap(CFStore2& cfs,
00141 ConvolutionFunction& cf,
00142 const VisBuffer2& vb, const Quantity& dPA,
00143 const Vector<Int>& dataChan2ImChanMap,
00144 const Vector<Int>& dataPol2ImPolMap,
00145 const Vector<Double>& pointingOffset);
00146
00147 Double runTimeG_p, runTimeDG_p, runTimeG1_p, runTimeG2_p, runTimeG3_p, runTimeG4_p, runTimeG5_p, runTimeG6_p, runTimeG7_p;
00148 Timer timer_p;
00149
00150
00151
00152
00153
00154
00155 protected:
00156 Vector<Double> uvwScale_p, offset_p, dphase_p;
00157 Vector<Int> chanMap_p, polMap_p;
00158 Matrix<Double> spwChanFreq_p, spwChanConjFreq_p;
00159 CFStore convFuncStore_p;
00160
00161 Vector<Int> inc_p;
00162 Int* __restrict__ incPtr_p;
00163 Vector<Int> cfMap_p, conjCFMap_p;
00164 VBRow2CFMapType vbRow2CFMap_p;
00165 VBRow2CFBMapType vbRow2CFBMap_p;
00166
00167
00168 void sgrid(Int& ndim,
00169 Double* __restrict__ pos,
00170 Int* __restrict__ loc,
00171 Int* __restrict__ off,
00172 Complex& phasor, const Int& irow,
00173 const Double* __restrict__ uvw,
00174 const Double& dphase, const Double& freq,
00175 const Double* __restrict__ scale,
00176 const Double* __restrict__ offset,
00177 const Float* __restrict__ sampling);
00178
00179 inline Bool onGrid (const Int& nx, const Int& ny,
00180 const Vector<Int>& __restrict__ loc,
00181 const Vector<Int>& __restrict__ support) __restrict__
00182 {
00183 return (((loc(0)-support[0]) >= 0 ) && ((loc(0)+support[0]) < nx) &&
00184 ((loc(1)-support[1]) >= 0 ) && ((loc(1)+support[1]) < ny));
00185 };
00186 inline Bool onGrid (const Int& nx, const Int& ny,
00187 const Int& loc0, const Int& loc1,
00188 const Int& support) __restrict__
00189 {
00190 return (((loc0-support) >= 0 ) && ((loc0+support) < nx) &&
00191 ((loc1-support) >= 0 ) && ((loc1+support) < ny));
00192 };
00193
00194
00195
00196
00197
00198
00199
00200
00201
00202
00203
00204
00205
00206
00207
00208
00209 inline void cacheAxisIncrements(const Int& n0, const Int& n1, const Int& n2, const Int& n3)
00210 {
00211
00212 inc_p.resize(4);
00213 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;
00214 Bool D;
00215 incPtr_p = inc_p.getStorage(D);
00216 }
00217 inline void cacheAxisIncrements(const Vector<Int>& n)
00218 {cacheAxisIncrements(n[0],n[1],n[2],n[3]);}
00219
00220 inline void cacheAxisIncrements(const Vector<Int>& n, Vector<Int>& inc)
00221 {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];}
00222
00223 inline void cacheAxisIncrements(const Int n[4], Int inc[4])
00224 {inc[0]=1; inc[1]=inc[0]*n[0]; inc[2]=inc[1]*n[1]; inc[3]=inc[2]*n[2];(void)n[3];}
00225
00226
00227
00228 inline void addTo4DArray(DComplex* __restrict__& store, Int* __restrict__& iPos,
00229 Complex& nvalue, Double& wt) __restrict__
00230 {addTo4DArray(store, iPos, incPtr_p, nvalue, wt);}
00231
00232 inline void addTo4DArray(Complex* __restrict__& store, Int* __restrict__& iPos,
00233 Complex& nvalue, Double& wt) __restrict__
00234 {addTo4DArray(store, iPos, incPtr_p, nvalue, wt);}
00235
00236
00237
00238 inline void addTo4DArray(DComplex* __restrict__& store, Int* __restrict__& iPos,
00239 Int* __restrict__ inc, Complex& nvalue, Double& wt) __restrict__
00240 {store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]] += (nvalue*Complex(wt));}
00241
00242 inline void addTo4DArray(Complex* __restrict__& store, Int* __restrict__& iPos,
00243 Int* __restrict__ inc, Complex& nvalue, Double& wt) __restrict__
00244 {store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]] += (nvalue*Complex(wt));}
00245
00246
00247 inline Complex getFrom4DArray(const Complex* __restrict__& store,
00248 const Int* __restrict__& iPos,
00249 const Vector<Int>& inc)
00250
00251 {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
00252
00253 inline Complex getFrom4DArray(const Complex* __restrict__& store,
00254 const Vector<Int> iPos, const Vector<Int>& inc)
00255
00256 {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
00257
00258 inline DComplex getFrom4DArray(const DComplex* __restrict__& store,
00259 const Int* __restrict__& iPos,
00260 const Vector<Int>& inc)
00261
00262 {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
00263
00264 inline DComplex getFrom4DArray(const DComplex* __restrict__& store,
00265 const Vector<Int> iPos, const Vector<Int>& inc)
00266
00267 {return store[iPos[0] + iPos[1]*inc[1] + iPos[2]*inc[2] +iPos[3]*inc[3]];};
00268
00269
00270
00271 inline Complex getFrom4DArray(const Complex* __restrict__& store, const Int* __restrict__& iPos)
00272
00273 {return getFrom4DArray(store, iPos, inc_p);}
00274
00275 inline DComplex getFrom4DArray(const DComplex* __restrict__& store, const Int* __restrict__& iPos)
00276
00277 {return getFrom4DArray(store, iPos, inc_p);}
00278
00279 };
00280 };
00281 };
00282 #endif //