Line data Source code
1 : //# GridFT.h: Definition for GridFT
2 : //# Copyright (C) 1996-2012
3 : //# Associated Universities, Inc. Washington DC, USA.
4 : //#
5 : //# This library is free software; you can redistribute it and/or modify it
6 : //# under the terms of the GNU Library General Public License as published by
7 : //# the Free Software Foundation; either version 2 of the License, or (at your
8 : //# option) any later version.
9 : //#
10 : //# This library is distributed in the hope that it will be useful, but WITHOUT
11 : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 : //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
13 : //# License for more details.
14 : //#
15 : //# You should have received a copy of the GNU Library General Public License
16 : //# along with this library; if not, write to the Free Software Foundation,
17 : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
18 : //#
19 : //# Correspondence concerning AIPS++ should be adressed as follows:
20 : //# Internet email: casa-feedback@nrao.edu.
21 : //# Postal address: AIPS++ Project Office
22 : //# National Radio Astronomy Observatory
23 : //# 520 Edgemont Road
24 : //# Charlottesville, VA 22903-2475 USA
25 : //#
26 : //#
27 : //# $Id$
28 :
29 : #ifndef SYNTHESIS_GRIDFT_H
30 : #define SYNTHESIS_GRIDFT_H
31 :
32 : #include <synthesis/TransformMachines/FTMachine.h>
33 : #include <casacore/casa/Arrays/Matrix.h>
34 : #include <casacore/scimath/Mathematics/FFTServer.h>
35 : #include <msvis/MSVis/VisBuffer.h>
36 : #include <casacore/images/Images/ImageInterface.h>
37 : #include <casacore/images/Images/ImageInterface.h>
38 : #include <casacore/casa/Containers/Block.h>
39 : #include <casacore/casa/Arrays/Array.h>
40 : #include <casacore/casa/Arrays/Vector.h>
41 : #include <casacore/casa/Arrays/Matrix.h>
42 : #include <casacore/scimath/Mathematics/ConvolveGridder.h>
43 : #include <casacore/lattices/Lattices/LatticeCache.h>
44 : #include <casacore/lattices/Lattices/ArrayLattice.h>
45 : //#include <synthesis/MeasurementComponents/SynthesisPeek.h>
46 :
47 :
48 : namespace casacore{
49 :
50 : class UVWMachine;
51 : }
52 :
53 : namespace casa { //# NAMESPACE CASA - BEGIN
54 :
55 : // <summary> An FTMachine for Gridded Fourier transforms </summary>
56 :
57 : // <use visibility=export>
58 :
59 : // <reviewed reviewer="" date="" tests="" demos="">
60 :
61 : // <prerequisite>
62 : // <li> <linkto class=FTMachine>FTMachine</linkto> module
63 : // <li> <linkto class=SkyEquation>SkyEquation</linkto> module
64 : // <li> <linkto class=VisBuffer>VisBuffer</linkto> module
65 : // </prerequisite>
66 : //
67 : // <etymology>
68 : // FTMachine is a Machine for Fourier Transforms. GridFT does
69 : // Grid-based Fourier transforms.
70 : // </etymology>
71 : //
72 : // <synopsis>
73 : // The <linkto class=SkyEquation>SkyEquation</linkto> needs to be able
74 : // to perform Fourier transforms on visibility data. GridFT
75 : // allows efficient Fourier Transform processing using a
76 : // <linkto class=VisBuffer>VisBuffer</linkto> which encapsulates
77 : // a chunk of visibility (typically all baselines for one time)
78 : // together with all the information needed for processing
79 : // (e.g. UVW coordinates).
80 : //
81 : // Gridding and degridding in GridFT are performed using a
82 : // novel sort-less algorithm. In this approach, the gridded plane is
83 : // divided into small patches, a cache of which is maintained in memory
84 : // using a general-purpose <linkto class=casacore::LatticeCache>LatticeCache</linkto> class. As the (time-sorted)
85 : // visibility data move around slowly in the Fourier plane, patches are
86 : // swapped in and out as necessary. Thus, optimally, one would keep at
87 : // least one patch per baseline.
88 : //
89 : // A grid cache is defined on construction. If the gridded uv plane is smaller
90 : // than this, it is kept entirely in memory and all gridding and
91 : // degridding is done entirely in memory. Otherwise a cache of tiles is
92 : // kept an paged in and out as necessary. Optimally the cache should be
93 : // big enough to hold all polarizations and frequencies for all
94 : // baselines. The paging rate will then be small. As the cache size is
95 : // reduced below this critical value, paging increases. The algorithm will
96 : // work for only one patch but it will be very slow!
97 : //
98 : // This scheme works well for arrays having a moderate number of
99 : // antennas since the saving in space goes as the ratio of
100 : // baselines to image size. For the ATCA, VLBA and WSRT, this ratio is
101 : // quite favorable. For the VLA, one requires images of greater than
102 : // about 200 pixels on a side to make it worthwhile.
103 : //
104 : // The FFT step is done plane by plane for images having less than
105 : // 1024 * 1024 pixels on each plane, and line by line otherwise.
106 : //
107 : // The gridding and degridding steps are implemented in Fortran
108 : // for speed. In gridding, the visibilities are added onto the
109 : // grid points in the neighborhood using a weighting function.
110 : // In degridding, the value is derived by a weight summ of the
111 : // same points, using the same weighting function.
112 : // </synopsis>
113 : //
114 : // <example>
115 : // See the example for <linkto class=SkyModel>SkyModel</linkto>.
116 : // </example>
117 : //
118 : // <motivation>
119 : // Define an interface to allow efficient processing of chunks of
120 : // visibility data
121 : // </motivation>
122 : //
123 : // <todo asof="97/10/01">
124 : // <ul> Deal with large VLA spectral line case
125 : // </todo>
126 :
127 : class GridFT : public FTMachine {
128 : public:
129 :
130 : // Constructor: cachesize is the size of the cache in words
131 : // (e.g. a few million is a good number), tilesize is the
132 : // size of the tile used in gridding (cannot be less than
133 : // 12, 16 works in most cases), and convType is the type of
134 : // gridding used (SF is prolate spheriodal wavefunction,
135 : // and BOX is plain box-car summation). mLocation is
136 : // the position to be used in some phase rotations. If
137 : // mTangent is specified then the uvw rotation is done for
138 : // that location iso the image center.
139 : // <group>
140 : GridFT();
141 : GridFT(casacore::Long cachesize, casacore::Int tilesize, casacore::String convType="SF",
142 : casacore::Float padding=1.0, casacore::Bool usezero=true, casacore::Bool useDoublePrec=false);
143 : GridFT(casacore::Long cachesize, casacore::Int tilesize, casacore::String convType,
144 : casacore::MPosition mLocation, casacore::Float padding=1.0, casacore::Bool usezero=true,
145 : casacore::Bool useDoublePrec=false);
146 : GridFT(casacore::Long cachesize, casacore::Int tilesize, casacore::String convType,
147 : casacore::MDirection mTangent, casacore::Float padding=1.0, casacore::Bool usezero=true,
148 : casacore::Bool useDoublePrec=false);
149 : GridFT(casacore::Long cachesize, casacore::Int tilesize, casacore::String convType,
150 : casacore::MPosition mLocation, casacore::MDirection mTangent, casacore::Float passing=1.0,
151 : casacore::Bool usezero=true, casacore::Bool useDoublePrec=false);
152 : // </group>
153 :
154 : // Construct from a casacore::Record containing the GridFT state
155 : GridFT(const casacore::RecordInterface& stateRec);
156 :
157 : // Copy constructor
158 : GridFT(const GridFT &other);
159 :
160 : // Assignment operator
161 : virtual GridFT &operator=(const GridFT &other);
162 :
163 : virtual ~GridFT();
164 :
165 : virtual FTMachine* cloneFTM();
166 :
167 : // Initialize transform to Visibility plane using the image
168 : // as a template. The image is loaded and Fourier transformed.
169 : virtual void initializeToVis(casacore::ImageInterface<casacore::Complex>& image,
170 : const VisBuffer& vb);
171 :
172 : // Finalize transform to Visibility plane: flushes the image
173 : // cache and shows statistics if it is being used.
174 : virtual void finalizeToVis();
175 :
176 : // Initialize transform to Sky plane: initializes the image
177 : virtual void initializeToSky(casacore::ImageInterface<casacore::Complex>& image, casacore::Matrix<casacore::Float>& weight,
178 : const VisBuffer& vb);
179 :
180 :
181 : // Finalize transform to Sky plane: flushes the image
182 : // cache and shows statistics if it is being used. DOES NOT
183 : // DO THE FINAL TRANSFORM!
184 : virtual void finalizeToSky();
185 :
186 :
187 : // Get actual coherence from grid by degridding
188 : virtual void get(VisBuffer& vb, casacore::Int row=-1);
189 :
190 :
191 : // Put coherence to grid by gridding.
192 : virtual void put(const VisBuffer& vb, casacore::Int row=-1, casacore::Bool dopsf=false,
193 : FTMachine::Type type=FTMachine::OBSERVED);
194 :
195 : // Make the entire image
196 : void makeImage(FTMachine::Type type,
197 : VisSet& vs,
198 : casacore::ImageInterface<casacore::Complex>& image,
199 : casacore::Matrix<casacore::Float>& weight);
200 :
201 : // Get the final image: do the Fourier transform and
202 : // grid-correct, then optionally normalize by the summed weights
203 : casacore::ImageInterface<casacore::Complex>& getImage(casacore::Matrix<casacore::Float>&, casacore::Bool normalize=true);
204 0 : virtual void normalizeImage(casacore::Lattice<casacore::Complex>& /*skyImage*/,
205 : const casacore::Matrix<casacore::Double>& /*sumOfWts*/,
206 : casacore::Lattice<casacore::Float>& /*sensitivityImage*/,
207 : casacore::Bool /*fftNorm*/)
208 0 : {throw(casacore::AipsError("GridFT::normalizeImage() called"));}
209 :
210 : // Get the final weights image
211 : void getWeightImage(casacore::ImageInterface<casacore::Float>&, casacore::Matrix<casacore::Float>&);
212 :
213 : // Save and restore the GridFT to and from a record
214 : virtual casacore::Bool toRecord(casacore::String& error, casacore::RecordInterface& outRec,
215 : casacore::Bool withImage=false, const casacore::String diskimage="");
216 : virtual casacore::Bool fromRecord(casacore::String& error, const casacore::RecordInterface& inRec);
217 :
218 : // Can this FTMachine be represented by Fourier convolutions?
219 0 : virtual casacore::Bool isFourier() {return true;}
220 :
221 0 : virtual void setNoPadding(casacore::Bool nopad){noPadding_p=nopad;};
222 : virtual void modifyConvFunc(const casacore::Vector<casacore::Double>& convFunc, casacore::Int convSupport, casacore::Int convSampling);
223 : virtual casacore::String name() const;
224 0 : virtual void setMiscInfo(const casacore::Int qualifier){(void)qualifier;};
225 0 : virtual void ComputeResiduals(VisBuffer&/*vb*/, casacore::Bool /*useCorrected*/) {};
226 :
227 : protected:
228 :
229 :
230 : // Padding in FFT
231 : casacore::Float padding_p;
232 :
233 : // Get the appropriate data pointer
234 : casacore::Array<casacore::Complex>* getDataPointer(const casacore::IPosition&, casacore::Bool);
235 :
236 : virtual void ok();
237 :
238 : virtual void init();
239 :
240 : //Prepare the grid for degridding
241 : virtual void prepGridForDegrid();
242 :
243 : // Is this record on Grid? check both ends. This assumes that the
244 : // ends bracket the middle
245 : casacore::Bool recordOnGrid(const VisBuffer& vb, casacore::Int rownr) const;
246 :
247 :
248 : // Image cache
249 : casacore::LatticeCache<casacore::Complex> * imageCache;
250 :
251 : // Sizes
252 : casacore::Long cachesize;
253 : casacore::Int tilesize;
254 :
255 : // Gridder
256 : casacore::ConvolveGridder<casacore::Double, casacore::Complex>* gridder;
257 :
258 : // Is this tiled?
259 : casacore::Bool isTiled;
260 :
261 : // casacore::Array lattice
262 : casacore::CountedPtr<casacore::Lattice<casacore::Complex> > arrayLattice;
263 :
264 : // Lattice. For non-tiled gridding, this will point to arrayLattice,
265 : // whereas for tiled gridding, this points to the image
266 : casacore::CountedPtr<casacore::Lattice<casacore::Complex> > lattice;
267 :
268 : casacore::String convType;
269 :
270 : casacore::Float maxAbsData;
271 :
272 : // Useful IPositions
273 : casacore::IPosition centerLoc, offsetLoc;
274 :
275 : // Image Scaling and offset
276 : casacore::Vector<casacore::Double> uvScale, uvOffset;
277 :
278 :
279 : casacore::Int priorCacheSize;
280 :
281 : // Grid/degrid zero spacing points?
282 :
283 : casacore::Bool usezero_p;
284 :
285 : //force no padding
286 : casacore::Bool noPadding_p;
287 :
288 : //Check if using put that avoids non-necessary reads
289 : casacore::Bool usePut2_p;
290 :
291 : //machine name
292 : casacore::String machineName_p;
293 :
294 : casacore::Double timemass_p, timegrid_p, timedegrid_p;
295 : casacore::Vector<casacore::Double> convFunc_p;
296 : casacore::Int convSampling_p, convSupport_p;
297 : // casa::async::SynthesisAsyncPeek *peek;
298 : };
299 :
300 : } //# NAMESPACE CASA - END
301 :
302 : #endif
|
