LSSTApplications  10.0-2-g4f67435,11.0.rc2+1,11.0.rc2+12,11.0.rc2+3,11.0.rc2+4,11.0.rc2+5,11.0.rc2+6,11.0.rc2+7,11.0.rc2+8
LSSTDataManagementBasePackage
offsetImage.cc
Go to the documentation of this file.
1 // -*- lsst-c++ -*-
2 
3 /*
4  * LSST Data Management System
5  * Copyright 2008, 2009, 2010 LSST Corporation.
6  *
7  * This product includes software developed by the
8  * LSST Project (http://www.lsst.org/).
9  *
10  * This program is free software: you can redistribute it and/or modify
11  * it under the terms of the GNU General Public License as published by
12  * the Free Software Foundation, either version 3 of the License, or
13  * (at your option) any later version.
14  *
15  * This program is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  * GNU General Public License for more details.
19  *
20  * You should have received a copy of the LSST License Statement and
21  * the GNU General Public License along with this program. If not,
22  * see <http://www.lsstcorp.org/LegalNotices/>.
23  */
24 
30 #include <iterator>
32 #include "lsst/afw/geom/Box.h"
33 #include "lsst/afw/geom/Extent.h"
36 
37 namespace afwImage = lsst::afw::image;
38 namespace afwGeom = lsst::afw::geom;
39 
40 namespace lsst {
41 namespace afw {
42 namespace math {
43 
54 template<typename ImageT>
55 typename ImageT::Ptr offsetImage(ImageT const& inImage,
56  float dx,
57  float dy,
58  std::string const& algorithmName,
59  unsigned int buffer
60  ) {
65  SeparableKernel::Ptr offsetKernel = makeWarpingKernel(algorithmName);
66 
67  typename ImageT::Ptr buffImage;
68  if (buffer > 0) {
69  // Paste input image into buffered image
70  afwGeom::Extent2I const &dims = inImage.getDimensions();
71  typename ImageT::Ptr buffered(new ImageT(dims.getX() + 2 * buffer, dims.getY() + 2 * buffer));
72  buffImage = buffered;
73  afwGeom::Box2I box(afwGeom::Point2I(buffer, buffer), dims);
74  typename ImageT::Ptr buffSmall(new ImageT(*buffImage, box, afwImage::LOCAL, false));
75  *buffSmall <<= inImage;
76  } else {
77  buffImage = boost::make_shared<ImageT>(inImage);
78  }
79 
80  if (offsetKernel->getWidth() > buffImage->getWidth() ||
81  offsetKernel->getHeight() > buffImage->getHeight()) {
82  throw LSST_EXCEPT(pexExcept::LengthError,
83  (boost::format("Image of size %dx%d is too small to offset using a %s kernel"
84  "(minimum %dx%d)") %
85  buffImage->getWidth() % buffImage->getHeight() % algorithmName %
86  offsetKernel->getWidth() % offsetKernel->getHeight()).str());
87  }
88 
89 // typename ImageT::Ptr convImage(new ImageT(buffImage, true)); // output image, a deep copy
90  typename ImageT::Ptr convImage(new ImageT(buffImage->getDimensions())); // Convolved image
91 
92  int dOrigX, dOrigY;
93  double fracX, fracY;
94  // If the offset in both axes is in (-1, 1) use it as is, and don't shift the origin
95  if (dx > -1 && dx < 1 && dy > -1 && dy < 1) {
96  dOrigX = 0;
97  dOrigY = 0;
98  fracX = dx;
99  fracY = dy;
100  } else {
101  dOrigX = static_cast<int>(std::floor(dx + 0.5));
102  dOrigY = static_cast<int>(std::floor(dy + 0.5));
103  fracX = dx - dOrigX;
104  fracY = dy - dOrigY;
105  }
106 
107  // We seem to have to pass -fracX, -fracY to setKernelParameters, for reasons RHL doesn't understand
108  double dKerX = -fracX;
109  double dKerY = -fracY;
110 
111  //
112  // If the shift is -ve, the generated shift kernel (e.g. Lanczos5) is quite asymmetric, with the
113  // largest coefficients to the left of centre. We therefore move the centre of calculated shift kernel
114  // one to the right to center up the largest coefficients
115  //
116  if (dKerX < 0) {
117  offsetKernel->setCtrX(offsetKernel->getCtrX() + 1);
118  }
119  if (dKerY < 0) {
120  offsetKernel->setCtrY(offsetKernel->getCtrY() + 1);
121  }
122 
123  offsetKernel->setKernelParameters(std::make_pair(dKerX, dKerY));
124 
125  convolve(*convImage, *buffImage, *offsetKernel, true, true);
126 
127  typename ImageT::Ptr outImage;
128  if (buffer > 0) {
129  afwGeom::Box2I box(afwGeom::Point2I(buffer, buffer), inImage.getDimensions());
130  typename ImageT::Ptr out(new ImageT(*convImage, box, afwImage::LOCAL, true));
131  outImage = out;
132  } else {
133  outImage = convImage;
134  }
135 
136  // adjust the origin; do this after convolution since convolution also sets XY0
137  outImage->setXY0(geom::Point2I(inImage.getX0() + dOrigX, inImage.getY0() + dOrigY));
138 
139  return outImage;
140 }
141 
142 /************************************************************************************************************/
143 //
144 // Explicit instantiations
145 //
147 #define INSTANTIATE(TYPE) \
148  template afwImage::Image<TYPE>::Ptr offsetImage(afwImage::Image<TYPE> const&, float, float, \
149  std::string const&, unsigned int); \
150  template afwImage::MaskedImage<TYPE>::Ptr offsetImage(afwImage::MaskedImage<TYPE> const&, float, float, \
151  std::string const&, unsigned int);
152 
153 INSTANTIATE(double)
154 INSTANTIATE(float)
155 INSTANTIATE(int)
157 
158 }}}
boost::shared_ptr< SeparableKernel > makeWarpingKernel(std::string name)
Return a warping kernel given its name.
Image utility functions.
An integer coordinate rectangle.
Definition: Box.h:53
table::Key< table::Array< Kernel::Pixel > > image
Definition: FixedKernel.cc:117
ImageT::Ptr offsetImage(ImageT const &image, float dx, float dy, std::string const &algorithmName="lanczos5", unsigned int buffer=0)
Return an image offset by (dx, dy) using the specified algorithm.
Definition: offsetImage.cc:55
#define LSST_EXCEPT(type,...)
Definition: Exception.h:46
#define INSTANTIATE(T)
Extent< int, N > floor(Extent< double, N > const &input)
void convolve(OutImageT &convolvedImage, InImageT const &inImage, KernelT const &kernel, ConvolutionControl const &convolutionControl=ConvolutionControl())
Convolve an Image or MaskedImage with a Kernel, setting pixels of an existing output image...
Support for warping an image to a new WCS.
A coordinate class intended to represent offsets and dimensions.
boost::shared_ptr< SeparableKernel > Ptr
Definition: Kernel.h:988