DM-10429: Write a version of the warper that uses SkyWcs and compare performance #228

r-owen · 2017-05-15T15:23:37Z

No description provided.

pschella

Looks good, except that the interpolation section is a bit long.

pschella · 2017-05-15T15:33:31Z

examples/timeWarpExposureUsingTransform.py

+
+#
+# LSST Data Management System
+# Copyright 2008, 2009, 2010 LSST Corporation.


Update copyright

pschella · 2017-05-15T15:35:03Z

examples/timeWarpExposureUsingTransform.py

+    int warpImage(
+        DestMaskedImageT &destMaskedImage,
+        SrcMaskedImageT const &srcMaskedImage,
+        SeparableKernel &warpingKernel, int const interpLength=0);


What is this C++ signature doing in Python code? Is this common?

Left in by accident

pschella · 2017-05-15T15:39:00Z

examples/timeWarpExposureUsingTransform.py

+def run():
+    if len(sys.argv) < 2:
+        srcExposure = afwImage.ExposureF(InputExposurePath)
+        if True:


Why this check? Just to have block scope? If so please indicate.

It intended to allow warping the whole image (slower) instead of a subregion. I added a global constant to this and timeWarpExposure.py to clarify this.

pschella · 2017-05-15T15:42:55Z

examples/timeWarpExposureUsingTransform.py

+    for interpLength in (0, 1, 5, 10):
+        for scaleFac in (1.2,):
+            destScale = srcScale / scaleFac
+            for offset_deg_arcsec in ((0.0, 0.0),):  # ((0.0, 0.0), (-35.0, 10.5)):


Uses different convention for variable naming.

pschella · 2017-05-15T15:49:17Z

include/lsst/afw/math/warpExposure.h

+ * @param[in,out] destImage  Destination image; all pixels are set
+ * @param[in] srcImage  Source image
+ * @param[in] srcFromDest  Transformation from destination to source parent pixel position
+ *    This can be computed as srcWcs.getInverse().of(destWcs)


I find this confusing, but that is probably just due to a lack of understanding on my part.

Why source parent instead of just source?
I would find something like srcWcs.transformationFrom(destWcs) more readable. But that might be incorrect.

Is also confusing in the test.

I'm afraid that's how it's done. Warping requires a transformation from destination pixel position to source. I could name it destToSrc but the same math is still required to derive it. Once our old image.Wcs class is replaced with geom.SkyWcs you can call the override of warpImage that takes the source and destination WCS directly, or call warpExposure, and not have to think about it at all.

pschella · 2017-05-15T16:11:46Z

src/math/warpExposure.cc

+        }
+        return 0;
+    }
+    std::shared_ptr<SeparableKernel> warpingKernelPtr = control.getWarpingKernel();


Why get it twice?

Or is it this call that can fail?

pschella · 2017-05-15T16:18:05Z

src/math/warpExposure.cc

+
+    int const destWidth = destImage.getWidth();
+    int const destHeight = destImage.getHeight();
+


remove (inconsistent) newline

pschella · 2017-05-15T16:40:07Z

tests/testWarpExposure.py

+            warpedSkyWcs = afwGeom.SkyWcs(swarpedMetadata)
+
+            # original image is source, warped image is destination
+            srcFromDest = originalSkyWcs.getInverse().of(warpedSkyWcs)


Again, this is confusing. See also the previous comments.

If we want it to be less confusing I fear some syntactic sugar may be required. I'm reluctant to fatten the API of SkyWcs until we are further along, preferring to document what's going on, especially since (as I say above) few users will ever need to do this. The overload being tested will not be the main way anyone warps things when they are using WCS to specify the warping. It will primary be used for simple transformations such as shifts and rotations, and as the implementation for the usual overrides.

pschella · 2017-05-15T16:42:02Z

tests/testWarpExposure.py

+
+            numGoodPix = afwMath.warpImage(afwWarpedMaskedImage, originalMaskedImage,
+                                           srcFromDest, warpingControl)
+            self.assertGreater(numGoodPix, 50)


To make sure that some pixels were warped.

pschella · 2017-05-15T18:29:13Z

src/math/warpExposure.cc

+
+    detail::WarpAtOnePoint<DestImageT, SrcImageT> warpAtOnePoint(srcImage, control, padValue);
+
+    if (interpLength > 0) {


The part in this branch can definitely benefit from some refactoring.
It looks ok, but there is way too much going on, I would suggest factoring some steps out into helper functions.

I agree, but this situation is temporary and I couldn't see an easy way to share code between the two blocks. Once SkyWcs and Transform are adopted, the similar code that handles XYTransform will disappear, leaving only the new code.

I wasn't thinking about sharing code between the blocks, but rather separating part of the one long block out a bit. It really is quite long. But there is no pressing need for it. Just feel like it would make it easier to grok.

I'll take a stab at it.

@pschella I pushed an additional commit with some refactoring of the interpolated warping code. It's not a big change, but it cleans up some code and avoids messing with tight loops where performance is a concern.

I want to withdraw my most recent change. It needlessly caches row information that can be computed on the fly. A clever implementation could compute the end row/column of each band and the fractional width on the fly (storing almost nothing) because all bands except the last are the same width. But that's more cleverness than I want to code.

pschella · 2017-05-16T01:02:06Z

src/math/warpExposure.cc

+            }  // for col
+            // move points from srcPosList to prevSrcPosList (we don't care about what ends up in srcPosList
+            // because it will be reallocated anyway)
+            std::swap(srcPosList, prevSrcPosList);


Generally prefer using std::swap to support user defined swap operations. In this case it is a std::vector so it doesn't matter, but we might change the type later :)

OK. I hate how that makes it hard to tell where swap is defined, but it does appear to be standard.

Added Endpoint.h and Transform.h to geom.h

Replace "if True" with a global variable to control whether `timeWarpExposure.py warps a subregion or the whole image

Add an overload of warpImage that uses a Transform<Point2Endpoint, Point2Endpoint> to describe the transform.

pschella reviewed May 15, 2017

View reviewed changes

r-owen force-pushed the tickets/DM-10429 branch 4 times, most recently from 6371990 to 94d193b Compare May 16, 2017 00:06

pschella reviewed May 16, 2017

View reviewed changes

r-owen added 4 commits May 17, 2017 10:16

Add missing imports to geom.h

950c598

Added Endpoint.h and Transform.h to geom.h

Add reference to makeCdMatrix to SkyWcs ctor

de3d685

Add global to control use subregion

d551213

Replace "if True" with a global variable to control whether `timeWarpExposure.py warps a subregion or the whole image

Add Transform-based overload of warpImage

05b99cd

Add an overload of warpImage that uses a Transform<Point2Endpoint, Point2Endpoint> to describe the transform.

r-owen force-pushed the tickets/DM-10429 branch from 3df10a0 to 05b99cd Compare May 17, 2017 17:40

r-owen merged commit 05b99cd into master May 17, 2017

ktlim deleted the tickets/DM-10429 branch August 25, 2018 06:43

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

DM-10429: Write a version of the warper that uses SkyWcs and compare performance #228

DM-10429: Write a version of the warper that uses SkyWcs and compare performance #228

r-owen commented May 15, 2017

pschella left a comment

pschella May 15, 2017

pschella May 15, 2017

r-owen May 15, 2017

pschella May 15, 2017

r-owen May 15, 2017

pschella May 15, 2017

pschella May 15, 2017

pschella May 15, 2017

pschella May 15, 2017

r-owen May 15, 2017

pschella May 15, 2017

pschella May 15, 2017

pschella May 15, 2017

pschella May 15, 2017

r-owen May 15, 2017 •

edited

pschella May 15, 2017

r-owen May 15, 2017

pschella May 15, 2017

r-owen May 15, 2017

pschella May 16, 2017

r-owen May 16, 2017

r-owen May 16, 2017

r-owen May 17, 2017 •

edited

pschella May 16, 2017

r-owen May 16, 2017


		int const destWidth = destImage.getWidth();
		int const destHeight = destImage.getHeight();


		detail::WarpAtOnePoint<DestImageT, SrcImageT> warpAtOnePoint(srcImage, control, padValue);

		if (interpLength > 0) {

DM-10429: Write a version of the warper that uses SkyWcs and compare performance #228

DM-10429: Write a version of the warper that uses SkyWcs and compare performance #228

Conversation

r-owen commented May 15, 2017

pschella left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

r-owen May 15, 2017 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

r-owen May 17, 2017 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

r-owen May 15, 2017 •

edited

r-owen May 17, 2017 •

edited