casacore
ImageRegrid.h
Go to the documentation of this file.
1 //# ImageRegrid.h: Regrid Images
2 //# Copyright (C) 1996,1997,1998,1999,2000,2001,2002,2003
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 addressed as follows:
20 //# Internet email: aips2-request@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 IMAGES_IMAGEREGRID_H
30 #define IMAGES_IMAGEREGRID_H
31 
32 #include <casacore/casa/aips.h>
33 #include <casacore/casa/Arrays/Matrix.h>
34 #include <casacore/casa/Arrays/Cube.h>
35 #include <casacore/measures/Measures/MDirection.h>
36 #include <casacore/measures/Measures/MFrequency.h>
37 #include <casacore/scimath/Mathematics/Interpolate2D.h>
38 #include <casacore/scimath/Mathematics/NumericTraits.h>
39 #include <set>
40 
41 namespace casacore { //# NAMESPACE CASACORE - BEGIN
42 
43 template<class T> class MaskedLattice;
44 template<class T> class ImageInterface;
45 template<class T> class Lattice;
46 template<class T> class LatticeIterator;
47 
48 class CoordinateSystem;
49 class DirectionCoordinate;
50 class Coordinate;
51 class ObsInfo;
52 class IPosition;
53 class Unit;
54 class ProgressMeter;
55 
56 // <summary>This regrids one image to match the coordinate system of another</summary>
57 
58 // <use visibility=export>
59 
60 // <reviewed reviewer="" date="yyyy/mm/dd" tests="" demos="">
61 // </reviewed>
62 
63 // <prerequisite>
64 // <li> <linkto class="ImageInterface">ImageInterface</linkto>
65 // <li> <linkto class="CoordinateSystem">CoordinateSystem</linkto>
66 // <li> <linkto class="Interpolate2D">Interpolate2D</linkto>
67 // <li> <linkto class="InterpolateArray1D">InterpolateArray1D</linkto>
68 // </prerequisite>
69 //
70 // <etymology>
71 // Regrids, or resamples, images.
72 // </etymology>
73 //
74 // <synopsis>
75 // This class enables you to regrid one image to the coordinate
76 // system of another. You can regrid any or all of the
77 // axes in the image. A range of interpolation schemes are available.
78 //
79 // It will cope with coordinate systems being in different orders
80 // (coordinate, world axes, pixel axes). The basic approach is to
81 // make a mapping from the input to the output coordinate systems,
82 // but the output CoordinateSystem order is preserved in the output
83 // image.
84 //
85 // Any DirectionCoordinate or LinearCoordinate holding exactly two axes
86 // is regridded in one pass with a 2-D interpolation scheme.
87 // All other axes are regridded in separate passes with a 1D interpolation
88 // scheme. This means that a LinearCoordinate holding say 3 axes
89 // where some of them are coupled will not be correctly regridded.
90 // StokesCoordinates cannot be regridded.
91 //
92 // Multiple passes are made through the data, and the output of
93 // each pass is the input of the next pass. The intermediate
94 // images are stored as TempImages which may be in memory or
95 // on disk, depending on their size.
96 //
97 // It can also simply insert this image into that one via
98 // an integer shift.
99 // </synopsis>
100 //
101 // <example>
102 //
103 // <srcblock>
104 // </srcblock>
105 // </example>
106 //
107 // <motivation>
108 // A common image analysis need is to regrid images, e.g. to compare
109 // images from different telescopes.
110 // </motivation>
111 //
112 // <thrown>
113 // <li> AipsError
114 // </thrown>
115 //
116 // <todo asof="1999/04/20">
117 // </todo>
118 
119 template <class T> class ImageRegrid
120 {
121 public:
122 
123  // Default constructor
125 
126  // copy constructor (copy semantics)
127  ImageRegrid(const ImageRegrid &other);
128 
129  // destructor
131 
132  // Assignment copy semantics)
134 
135  // Regrid inImage onto the grid specified by outImage.
136  // If outImage has a writable mask, it will be updated in that
137  // output pixels at which the regridding failed will be masked bad (False)
138  // and the pixel value set to zero. Otherwise the output mask is not changed.
139  // Specify which pixel axes of outImage are to be
140  // regridded. The coordinate and axis order of outImage
141  // is preserved, regardless of where the relevant coordinates
142  // are in inImage.
143  //
144  // decimate only applies when replicate=False. it is
145  // the coordinate grid computation decimation FACTOR
146  // (e.g. nCoordGrid ~ nIn / decimate). 0 means no decimation
147  // (slowest and most accurate)
148  void regrid(ImageInterface<T>& outImage,
149  typename Interpolate2D::Method method,
150  const IPosition& whichOutPixelAxes,
151  const ImageInterface<T>& inImage,
152  Bool replicate=False, uInt decimate=0,
153  Bool showProgress=False, Bool forceRegrid=False,
154  Bool verbose=False);
155 
156 // Get and set the 2-D coordinate grid. After a call to function <src>regrid</src>
157 // in which coupled 2D coordinate (presently only DirectionCoordinate) is
158 // regridded, this coordinate grid will be available. It can be reused
159 // via the <src>set2DCoordinateGrid</src> function for another like plane
160 // (e.g. if you choose to regrid planes of a cube separately). When you provide
161 // the coordinate grid, it will no longer (for that 2D coordinate only) be
162 // computed internally, which may save a lot of time. Ordinarily, if you
163 // regridded many planes of a cube in one call to regrid, the coordinate grid
164 // is cached for you. To trigger successive calls to regrid to go back to
165 // internal computation, set zero length Cube and Matrix. <src>gridMask</src>
166 // is True for successfull coordinate conversions, and False otherwise.
167 // <group>
168  void get2DCoordinateGrid (Cube<Double>& grid, Matrix<Bool>& gridMask) const;
169  void set2DCoordinateGrid (const Cube<Double>& grid, const Matrix<Bool>& gridMask, Bool notify=False);
170 // </group>
171 //
172  // Inserts inImage into outImage. The alignment is done by
173  // placing the blc of inImage at the specified
174  // absolute pixel of the outImage (outPixelLocation). If
175  // the outPixelLocation vector is of zero length, then the images
176  // are aligned by their reference pixels. Only integral shifts are done
177  // in the aligment process. If outImage has a mask, it will be updated.
178  // Returns False if no overlap of images, in which case the
179  // output is not updated.
181  const Vector<Double>& outPixelLocation,
182  const ImageInterface<T>& inImage);
183 
184  // Print out useful debugging information (level 0 is none,
185  // 1 is some, 2 is too much)
186  void showDebugInfo(Int level=0) {itsShowLevel = level;};
187 
188  // Enable/disable Measures Reference conversions
190 
191  // Helper function. We are regridding from cSysFrom to cSysTo for the
192  // specified pixel axes of cSyFrom. This function returns a CoordinateSystem which,
193  // for the pixel axes being regridded, copies the coordinates from cSysTo
194  // (if coordinate type present in cSysTo) or cSysFrom (coordinate
195  // type not present in cSysTo).
196  // For the axes not being regridded, it copies the coordinates from
197  // cSysFrom. This helps you build the cSys for function regrid.
198  // The ObsInfo from cSysFrom is copied to the output CoordinateSystem.
199  // If inShape has one or more elements it represenents the size of the
200  // image to be regridded. It this must have the same number of elements
201  // as the number of pixel axes in <src>cSysFrom</src>. If any of the values
202  // are unity (ie the axes are degenerate), and the corresponding axis in <src>csysFrom</src> is the only
203  // axis in its corresponding coordinate, this coordinate will not be replaced
204  // even if the axis is specified in <src>axes</src>.
205  // Upon return, <src>coordsToBeRegridded</src> will contain a list of the coordinates that will
206  // be regridded.
208  LogIO& os,
209  std::set<Coordinate::Type>& coordsToBeRegridded,
210  const CoordinateSystem& cSysTo,
211  const CoordinateSystem& cSysFrom,
212  const IPosition& axes,
213  const IPosition& inShape=IPosition(),
214  Bool giveStokesWarning=True
215  );
216 
217  private:
218 
221 //
224 //
228 //
229  // Check shape and axes. Exception if no good. If pixelAxes
230  // of length 0, set to all axes according to shape
231  void _checkAxes(IPosition& outPixelAxes,
232  const IPosition& inShape,
233  const IPosition& outShape,
234  const Vector<Int>& pixelAxisMap,
235  const CoordinateSystem& outCoords,
236  Bool verbose);
237 
238  // Find maps between coordinate systems
239  void findMaps (uInt nDim,
240  Vector<Int>& pixelAxisMap1,
241  Vector<Int>& pixelAxisMap2,
242  const CoordinateSystem& inCoords,
243  const CoordinateSystem& outCoords) const;
244 
245  // Find scale factor to conserve flux
246  Double findScaleFactor(const Unit& units,
247  const CoordinateSystem& inCoords,
248  const CoordinateSystem& outCoords,
249  Int inCoordinate, Int outCoordinate,
250  LogIO& os) const;
251 
252  // Regrid one Coordinate
253  void _regridOneCoordinate (LogIO& os, IPosition& outShape2,
254  Vector<Bool>& doneOutPixelAxes,
255  MaskedLattice<T>* &finalOutPtr,
256  MaskedLattice<T>* &inPtr,
257  MaskedLattice<T>* &outPtr,
258  CoordinateSystem& outCoords,
259  const CoordinateSystem& inCoords,
260  Int outPixelAxis,
261  const ImageInterface<T>& inImage,
262  const IPosition& outShape,
263  Bool replicate, uInt decimate,
264  Bool outIsMasked, Bool showProgress,
265  Bool forceRegrid,
266  typename Interpolate2D::Method method,
267  Bool verbose);
268 
269  // Regrid DirectionCoordinate or 2-axis LinearCoordinate
271  const MaskedLattice<T>& inLattice,
272  const Unit& imageUnit,
273  const CoordinateSystem& inCoords,
274  const CoordinateSystem& outCoords,
275  Int inCoordinate, Int outCoordinate,
276  const Vector<Int> inPixelAxes,
277  const Vector<Int> outPixelAxes,
278  const Vector<Int> pixelAxisMap1,
279  const Vector<Int> pixelAxisMap2,
280  typename Interpolate2D::Method method,
281  Bool replicate, uInt decimate,
282  Bool showProgress);
283 
284  // Make regridding coordinate grid for this cursor.
285  void make2DCoordinateGrid (LogIO& os, Bool& allFail, Bool&missedIt,
286  Double& minInX, Double& minInY,
287  Double& maxInX, Double& maxInY,
288  Cube<Double>& in2DPos,
289  Matrix<Bool>& succeed,
290  const CoordinateSystem& inCoords,
291  const CoordinateSystem& outCoords,
292  Int inCoordinate, Int outCoordinate,
293  uInt xInAxis, uInt yInAxis,
294  uInt xOutAxis, uInt yOutAxis,
295  const IPosition& inPixelAxes,
296  const IPosition& outPixelAxes,
297  const IPosition& inShape,
298  const IPosition& outPos,
299  const IPosition& cursorShape,
300  uInt decimate=0);
301 
302  // Make replication coordinate grid for this cursor
304  Double& minInX, Double& minInY,
305  Double& maxInX, Double& maxInY,
306  const Vector<Double>& pixelScale,
307  uInt xInAxis, uInt yInAxis,
308  uInt xOutAxis, uInt yOutAxis,
309  uInt xInCorrAxis, uInt yInCorrAxis,
310  uInt xOutCorrAxis, uInt yOutCorrAxis,
311  const IPosition& outPos, const IPosition& cursorShape);
312 
313  // Make regridding coordinate grid for this axis
315  Vector<Bool>& failed,
316  Bool& allFailed,
317  Bool& allGood,
318  const Coordinate& inCoord,
319  const Coordinate& outCoord,
320  Int inAxisInCoordinate,
321  Int outAxisInCoordinate,
322  MFrequency::Convert& machine,
323  Bool useMachine);
324 
325 
326  // Make replication coordinate grid for this axis
328  typename NumericTraits<T>::BaseType pixelScale) const;
329 
330  // Regrid 1 axis
331  void regrid1D (MaskedLattice<T>& outLattice,
332  const MaskedLattice<T>& inLattice,
333  const Coordinate& inCoord,
334  const Coordinate& outCoord,
335  const Vector<Int>& inPixelAxes,
336  const Vector<Int>& outPixelAxes,
337  Int inAxisInCoordinate,
338  Int outAxisInCoordinate,
339  const Vector<Int> pixelAxisMap,
340  typename Interpolate2D::Method method,
341  MFrequency::Convert& machine,
342  Bool replicate,
343  Bool useMachine, Bool showProgress);
344 
345 //
346  void regrid2DMatrix(Lattice<T>& outCursor,
347  LatticeIterator<Bool>*& outMaskIterPtr,
348  const Interpolate2D& interp,
349  ProgressMeter*& pProgress,
350  Double& iPix,
351  uInt nDim,
352  uInt xInAxis, uInt yInAxis,
353  uInt xOutAxis, uInt yOutAxis,
354  Double scale,
355  Bool inIsMasked, Bool outIsMasked,
356  const IPosition& outPos,
357  const IPosition& outCursorShape,
358  const IPosition& inChunkShape,
359  const IPosition& inChunkBlc,
360  const IPosition& pixelAxisMap2,
361  Array<T>& inDataChunk,
362  Array<Bool>*& inMaskChunkPtr,
363  const Cube<Double>& pix2DPos,
364  const Matrix<Bool>& succeed);
365 
366  void findXYExtent (Bool& missedIt, Bool& allFailed,
367  Double& minInX, Double& minInY,
368  Double& maxInX, Double& maxInY,
369  Cube<Double>& in2DPos,
370  const Matrix<Bool>& succeed,
371  uInt xInAxis, uInt yInAxis,
372  uInt xOutAxis, uInt yOutAxis,
373  const IPosition& outPos,
374  const IPosition& outCursorShape,
375  const IPosition& inShape);
376 //
377  Bool minmax(Double &minX, Double &maxX, Double& minY, Double& maxY,
378  const Array<Double> &xData,
379  const Array<Double> &yData,
380  const Array<Bool>& mask);
381 };
382 
383 //# Declare extern templates for often used types.
384  extern template class ImageRegrid<Float>;
385 
386 } //# NAMESPACE CASACORE - END
387 
388 #ifndef CASACORE_NO_AUTO_TEMPLATES
389 #include <casacore/images/Images/ImageRegrid.tcc>
390 #endif //# CASACORE_NO_AUTO_TEMPLATES
391 #endif
392 
simple 1-D array
Definition: Block.h:200
ImageRegrid()
Default constructor.
ImageRegrid(const ImageRegrid &other)
copy constructor (copy semantics)
Double findScaleFactor(const Unit &units, const CoordinateSystem &inCoords, const CoordinateSystem &outCoords, Int inCoordinate, Int outCoordinate, LogIO &os) const
Find scale factor to conserve flux
void get2DCoordinateGrid(Cube< Double > &grid, Matrix< Bool > &gridMask) const
Get and set the 2-D coordinate grid.
Matrix< Bool > itsUser2DCoordinateGridMask
Definition: ImageRegrid.h:226
Cube< Double > itsUser2DCoordinateGrid
Definition: ImageRegrid.h:225
void make2DCoordinateGrid(LogIO &os, Bool &allFail, Bool &missedIt, Double &minInX, Double &minInY, Double &maxInX, Double &maxInY, Cube< Double > &in2DPos, Matrix< Bool > &succeed, const CoordinateSystem &inCoords, const CoordinateSystem &outCoords, Int inCoordinate, Int outCoordinate, uInt xInAxis, uInt yInAxis, uInt xOutAxis, uInt yOutAxis, const IPosition &inPixelAxes, const IPosition &outPixelAxes, const IPosition &inShape, const IPosition &outPos, const IPosition &cursorShape, uInt decimate=0)
Make regridding coordinate grid for this cursor.
void regridTwoAxisCoordinate(LogIO &os, MaskedLattice< T > &outLattice, const MaskedLattice< T > &inLattice, const Unit &imageUnit, const CoordinateSystem &inCoords, const CoordinateSystem &outCoords, Int inCoordinate, Int outCoordinate, const Vector< Int > inPixelAxes, const Vector< Int > outPixelAxes, const Vector< Int > pixelAxisMap1, const Vector< Int > pixelAxisMap2, typename Interpolate2D::Method method, Bool replicate, uInt decimate, Bool showProgress)
Regrid DirectionCoordinate or 2-axis LinearCoordinate.
void regrid2DMatrix(Lattice< T > &outCursor, LatticeIterator< Bool > *&outMaskIterPtr, const Interpolate2D &interp, ProgressMeter *&pProgress, Double &iPix, uInt nDim, uInt xInAxis, uInt yInAxis, uInt xOutAxis, uInt yOutAxis, Double scale, Bool inIsMasked, Bool outIsMasked, const IPosition &outPos, const IPosition &outCursorShape, const IPosition &inChunkShape, const IPosition &inChunkBlc, const IPosition &pixelAxisMap2, Array< T > &inDataChunk, Array< Bool > *&inMaskChunkPtr, const Cube< Double > &pix2DPos, const Matrix< Bool > &succeed)
ImageRegrid< T > & operator=(const ImageRegrid &other)
Assignment copy semantics)
void make1DCoordinateGrid(Block< typename NumericTraits< T >::BaseType > &xOut, Vector< Bool > &failed, Bool &allFailed, Bool &allGood, const Coordinate &inCoord, const Coordinate &outCoord, Int inAxisInCoordinate, Int outAxisInCoordinate, MFrequency::Convert &machine, Bool useMachine)
Make regridding coordinate grid for this axis.
void regrid(ImageInterface< T > &outImage, typename Interpolate2D::Method method, const IPosition &whichOutPixelAxes, const ImageInterface< T > &inImage, Bool replicate=False, uInt decimate=0, Bool showProgress=False, Bool forceRegrid=False, Bool verbose=False)
Regrid inImage onto the grid specified by outImage.
void set2DCoordinateGrid(const Cube< Double > &grid, const Matrix< Bool > &gridMask, Bool notify=False)
Matrix< Bool > its2DCoordinateGridMask
Definition: ImageRegrid.h:223
void make2DCoordinateGrid(Cube< Double > &in2DPos, Double &minInX, Double &minInY, Double &maxInX, Double &maxInY, const Vector< Double > &pixelScale, uInt xInAxis, uInt yInAxis, uInt xOutAxis, uInt yOutAxis, uInt xInCorrAxis, uInt yInCorrAxis, uInt xOutCorrAxis, uInt yOutCorrAxis, const IPosition &outPos, const IPosition &cursorShape)
Make replication coordinate grid for this cursor.
void make1DCoordinateGrid(Block< typename NumericTraits< T >::BaseType > &xOut, typename NumericTraits< T >::BaseType pixelScale) const
Make replication coordinate grid for this axis.
void _checkAxes(IPosition &outPixelAxes, const IPosition &inShape, const IPosition &outShape, const Vector< Int > &pixelAxisMap, const CoordinateSystem &outCoords, Bool verbose)
Check shape and axes.
Bool minmax(Double &minX, Double &maxX, Double &minY, Double &maxY, const Array< Double > &xData, const Array< Double > &yData, const Array< Bool > &mask)
void findMaps(uInt nDim, Vector< Int > &pixelAxisMap1, Vector< Int > &pixelAxisMap2, const CoordinateSystem &inCoords, const CoordinateSystem &outCoords) const
Find maps between coordinate systems.
void findXYExtent(Bool &missedIt, Bool &allFailed, Double &minInX, Double &minInY, Double &maxInX, Double &maxInY, Cube< Double > &in2DPos, const Matrix< Bool > &succeed, uInt xInAxis, uInt yInAxis, uInt xOutAxis, uInt yOutAxis, const IPosition &outPos, const IPosition &outCursorShape, const IPosition &inShape)
void regrid1D(MaskedLattice< T > &outLattice, const MaskedLattice< T > &inLattice, const Coordinate &inCoord, const Coordinate &outCoord, const Vector< Int > &inPixelAxes, const Vector< Int > &outPixelAxes, Int inAxisInCoordinate, Int outAxisInCoordinate, const Vector< Int > pixelAxisMap, typename Interpolate2D::Method method, MFrequency::Convert &machine, Bool replicate, Bool useMachine, Bool showProgress)
Regrid 1 axis.
Cube< Double > its2DCoordinateGrid
Definition: ImageRegrid.h:222
void _regridOneCoordinate(LogIO &os, IPosition &outShape2, Vector< Bool > &doneOutPixelAxes, MaskedLattice< T > *&finalOutPtr, MaskedLattice< T > *&inPtr, MaskedLattice< T > *&outPtr, CoordinateSystem &outCoords, const CoordinateSystem &inCoords, Int outPixelAxis, const ImageInterface< T > &inImage, const IPosition &outShape, Bool replicate, uInt decimate, Bool outIsMasked, Bool showProgress, Bool forceRegrid, typename Interpolate2D::Method method, Bool verbose)
Regrid one Coordinate.
void disableReferenceConversions(Bool disable=True)
Enable/disable Measures Reference conversions.
Definition: ImageRegrid.h:189
~ImageRegrid()
destructor
void showDebugInfo(Int level=0)
Print out useful debugging information (level 0 is none, 1 is some, 2 is too much)
Definition: ImageRegrid.h:186
Bool insert(ImageInterface< T > &outImage, const Vector< Double > &outPixelLocation, const ImageInterface< T > &inImage)
Inserts inImage into outImage.
static CoordinateSystem makeCoordinateSystem(LogIO &os, std::set< Coordinate::Type > &coordsToBeRegridded, const CoordinateSystem &cSysTo, const CoordinateSystem &cSysFrom, const IPosition &axes, const IPosition &inShape=IPosition(), Bool giveStokesWarning=True)
Helper function.
A read/write lattice iterator.
Char BaseType
Numeric type.
this file contains all the compiler specific defines
Definition: mainpage.dox:28
const Bool False
Definition: aipstype.h:44
unsigned int uInt
Definition: aipstype.h:51
LatticeExprNode mask(const LatticeExprNode &expr)
This function returns the mask of the given expression.
int Int
Definition: aipstype.h:50
bool Bool
Define the standard types used by Casacore.
Definition: aipstype.h:42
const Bool True
Definition: aipstype.h:43
double Double
Definition: aipstype.h:55