mxlib
c++ tools for analyzing astronomical data and other tasks by Jared R. Males. [git repo]
Loading...
Searching...
No Matches
imageMasks_test.cpp
Go to the documentation of this file.
1/** \file imageMasks_test.cpp
2 * \brief Tests image-mask generation.
3 */
4#include "../../catch2/catch.hpp"
5
6#include <algorithm>
7#include <vector>
8#include <Eigen/Dense>
9
10#define MX_NO_ERROR_REPORTS
11
12#include "../../../include/improc/imageMasks.hpp"
13#include "../../../include/improc/eigenImage.hpp"
14
15/** \cond */
16/// Emit the binary-mask rotation specialization for LCOV accounting.
17template void mx::improc::rotateMask<mx::improc::eigenImage<double>, mx::improc::cubicConvolTransform<double>>(
18 mx::improc::eigenImage<double> &, mx::improc::eigenImage<double> &, double );
19/** \endcond */
20
21/** Masking wedges in an image
22 *
23 * Verify wedge masking, including that all pixels are masked for continuous rotations of the wedge
24 *
25 */
26/**
27 * \ingroup imageMasks_unit_tests
28 */
29TEST_CASE( "Masking wedges in an image", "[improc::imageMasks::maskWedge]" )
30{
31 GIVEN( "a single wedge" )
32 {
33 WHEN( "geometric center, 0-90 degrees" )
34 {
35 mx::improc::eigenImage<double> im;
36 im.resize( 1024, 1024 );
37 im.setZero();
38
39 double xcen = 0.5 * ( im.rows() - 1 );
40 double ycen = 0.5 * ( im.cols() - 1 );
41
42 mx::improc::maskWedge( im, xcen, ycen, 45.0, 45.0, 1 );
43
44 REQUIRE( im.sum() == 512 * 512 );
45 }
46 WHEN( "geometric center, 90-180 degrees" )
47 {
48 mx::improc::eigenImage<double> im;
49 im.resize( 1024, 1024 );
50 im.setZero();
51
52 double xcen = 0.5 * ( im.rows() - 1 );
53 double ycen = 0.5 * ( im.cols() - 1 );
54
55 mx::improc::maskWedge( im, xcen, ycen, 135.0, 45.0, 1 );
56
57 REQUIRE( im.sum() == 512 * 512 );
58 }
59 WHEN( "geometric center, 180-270 degrees" )
60 {
61 mx::improc::eigenImage<double> im;
62 im.resize( 1024, 1024 );
63 im.setZero();
64
65 double xcen = 0.5 * ( im.rows() - 1 );
66 double ycen = 0.5 * ( im.cols() - 1 );
67
68 mx::improc::maskWedge( im, xcen, ycen, 225.0, 45.0, 1 );
69
70 REQUIRE( im.sum() == 512 * 512 );
71 }
72 WHEN( "geometric center, 270-360 degrees" )
73 {
74 mx::improc::eigenImage<double> im;
75 im.resize( 1024, 1024 );
76 im.setZero();
77
78 double xcen = 0.5 * ( im.rows() - 1 );
79 double ycen = 0.5 * ( im.cols() - 1 );
80
81 mx::improc::maskWedge( im, xcen, ycen, 315.0, 45.0, 1 );
82
83 REQUIRE( im.sum() == 512 * 512 );
84 }
85 WHEN( "geometric center, 45-135 degrees" )
86 {
87 mx::improc::eigenImage<double> im;
88 im.resize( 1024, 1024 );
89 im.setZero();
90
91 double xcen = 0.5 * ( im.rows() - 1 );
92 double ycen = 0.5 * ( im.cols() - 1 );
93
94 mx::improc::maskWedge( im, xcen, ycen, 90.0, 45.0, 1 );
95
96 REQUIRE( im.sum() == 512 * 512 );
97 }
98 WHEN( "geometric center, 135-225 degrees" )
99 {
100 mx::improc::eigenImage<double> im;
101 im.resize( 1024, 1024 );
102 im.setZero();
103
104 double xcen = 0.5 * ( im.rows() - 1 );
105 double ycen = 0.5 * ( im.cols() - 1 );
106
107 mx::improc::maskWedge( im, xcen, ycen, 180.0, 45.0, 1 );
108
109 REQUIRE( im.sum() == 512 * 512 );
110 }
111 WHEN( "geometric center, 225-315 degrees" )
112 {
113 mx::improc::eigenImage<double> im;
114 im.resize( 1024, 1024 );
115 im.setZero();
116
117 double xcen = 0.5 * ( im.rows() - 1 );
118 double ycen = 0.5 * ( im.cols() - 1 );
119
120 mx::improc::maskWedge( im, xcen, ycen, 270.0, 45.0, 1 );
121
122 REQUIRE( im.sum() == 512 * 512 );
123 }
124 WHEN( "geometric center, 315-45 degrees" )
125 {
126 mx::improc::eigenImage<double> im;
127 im.resize( 1024, 1024 );
128 im.setZero();
129
130 double xcen = 0.5 * ( im.rows() - 1 );
131 double ycen = 0.5 * ( im.cols() - 1 );
132
133 mx::improc::maskWedge( im, xcen, ycen, 0.0, 45.0, 1 );
134
135 REQUIRE( im.sum() == 512 * 512 );
136 }
137 WHEN( "geometric center, 3 wedges of 120 degrees" )
138 {
139 mx::improc::eigenImage<double> im60;
140 im60.resize( 1024, 1024 );
141 im60.setZero();
142
143 double xcen = 0.5 * ( im60.rows() - 1 );
144 double ycen = 0.5 * ( im60.cols() - 1 );
145
146 mx::improc::maskWedge( im60, xcen, ycen, 60.0, 60.0, 1 );
147
148 mx::improc::eigenImage<double> im180;
149 im180.resize( 1024, 1024 );
150 im180.setZero();
151
152 xcen = 0.5 * ( im180.rows() - 1 );
153 ycen = 0.5 * ( im180.cols() - 1 );
154
155 mx::improc::maskWedge( im180, xcen, ycen, 180.0, 60.0, 1 );
156
157 mx::improc::eigenImage<double> im300;
158 im300.resize( 1024, 1024 );
159 im300.setZero();
160
161 xcen = 0.5 * ( im300.rows() - 1 );
162 ycen = 0.5 * ( im300.cols() - 1 );
163
164 mx::improc::maskWedge( im300, xcen, ycen, 300.0, 60.0, 1 );
165
166 REQUIRE( im60.sum() + im180.sum() + im300.sum() == 1024 * 1024 );
167 }
168 }
169}
170
171/** \brief Verifies mx::improc::rotateMask resizes, rotates, and thresholds binary masks.
172 *
173 * \ingroup imageMasks_unit_tests
174 */
175TEST_CASE( "rotateMask preserves binary masks through interpolation", "[improc::rotateMask]" )
176{
177 mx::improc::eigenImage<double> mask( 9, 9 );
178 mask.setZero();
179 mask( 3, 4 ) = 0.6;
180 mask( 4, 3 ) = 0.49;
181
182 mx::improc::eigenImage<double> rotated( 1, 1 );
183 mx::improc::rotateMask( rotated, mask, 0.0 );
184 REQUIRE( rotated.rows() == mask.rows() );
185 REQUIRE( rotated.cols() == mask.cols() );
186 REQUIRE( rotated( 3, 4 ) == 1.0 );
187 REQUIRE( rotated( 4, 3 ) == 0.0 );
188
189 mx::improc::rotateMask( rotated, mask, mx::math::half_pi<double>() );
190 REQUIRE( rotated( 4, 3 ) == 1.0 );
191 REQUIRE( rotated( 3, 4 ) == 0.0 );
192}
193
194/** \brief Verifies radius and degree-angle images for a centered rectangular grid.
195 *
196 * \ingroup imageMasks_unit_tests
197 */
198TEST_CASE( "radAngImage produces HCI radius and angle coordinates", "[improc::imageMasks::radAngImage]" )
199{
200 mx::improc::eigenImage<double> radius( 3, 5 );
201 mx::improc::eigenImage<double> angle;
202
203 mx::improc::radAngImage<mx::math::degreesT<double>>( radius, angle, 1.0, 2.0 );
204
205 REQUIRE( angle.rows() == 3 );
206 REQUIRE( angle.cols() == 5 );
207 REQUIRE( radius( 1, 2 ) == Approx( 0.0 ) );
208 REQUIRE( angle( 1, 2 ) == Approx( 0.0 ) );
209 REQUIRE( radius( 2, 2 ) == Approx( 1.0 ) );
210 REQUIRE( angle( 2, 2 ) == Approx( 0.0 ) );
211 REQUIRE( radius( 1, 3 ) == Approx( 1.0 ) );
212 REQUIRE( angle( 1, 3 ) == Approx( 90.0 ) );
213 REQUIRE( radius( 0, 2 ) == Approx( 1.0 ) );
214 REQUIRE( angle( 0, 2 ) == Approx( 180.0 ) );
215 REQUIRE( radius( 1, 1 ) == Approx( 1.0 ) );
216 REQUIRE( angle( 1, 1 ) == Approx( 270.0 ) );
217
218 mx::improc::radAngImage<mx::math::degreesT<double>>( radius, angle, 1.0, 2.0, 2.0 );
219 REQUIRE( radius( 2, 2 ) == Approx( 2.0 ) );
220}
221
222/** \brief Verifies HCI annulus indices for radial, angular, image-boundary, and mask constraints.
223 *
224 * \ingroup imageMasks_unit_tests
225 */
226TEST_CASE( "annulusIndices selects bounded and masked HCI regions", "[improc::imageMasks::annulusIndices]" )
227{
228 using angleT = mx::math::degreesT<double>;
229
230 mx::improc::eigenImage<double> radius( 5, 5 );
231 mx::improc::eigenImage<double> angle;
232 mx::improc::radAngImage<angleT>( radius, angle, 2.0, 2.0 );
233
234 const std::vector<size_t> full =
235 mx::improc::annulusIndices<angleT>( radius, angle, 2.0, 2.0, 1.0, 2.0, 0.0, 360.0 );
236 REQUIRE( full == std::vector<size_t>{ 6, 7, 8, 11, 13, 16, 17, 18 } );
237
238 mx::improc::eigenImage<double> mask( 5, 5 );
239 mask.setOnes();
240 mask( 1, 1 ) = 0;
241 const std::vector<size_t> masked =
242 mx::improc::annulusIndices<angleT>( radius, angle, 2.0, 2.0, 1.0, 2.0, 0.0, 360.0, &mask );
243 REQUIRE( masked == std::vector<size_t>{ 7, 8, 11, 13, 16, 17, 18 } );
244
245 const std::vector<size_t> clipped =
246 mx::improc::annulusIndices<angleT>( radius, angle, 2.0, 2.0, 0.0, 10.0, -360.0, 360.0 );
247 REQUIRE( clipped.size() == 25 );
248
249 const std::vector<size_t> half =
250 mx::improc::annulusIndices<angleT>( radius, angle, 2.0, 2.0, 0.0, 10.0, 0.0, 180.0 );
251 REQUIRE_FALSE( half.empty() );
252 REQUIRE( half.size() < clipped.size() );
253
254 const std::vector<size_t> wrapped =
255 mx::improc::annulusIndices<angleT>( radius, angle, 2.0, 2.0, 0.0, 10.0, 300.0, 60.0 );
256 REQUIRE_FALSE( wrapped.empty() );
257 REQUIRE( wrapped.size() < clipped.size() );
258}
259
260/** \brief Verifies annulusCoords and annulusIndices forward radial pixel buffers to their shared worker.
261 *
262 * \ingroup imageMasks_unit_tests
263 */
264TEST_CASE( "Annulus coordinate wrappers apply radial pixel buffers", "[improc::imageMasks][annulusCoords][pixbuf]" )
265{
266 using angleT = mx::math::degreesT<double>;
267
268 mx::improc::eigenImage<double> radius( 7, 7 );
269 mx::improc::eigenImage<double> angle;
270 mx::improc::radAngImage<angleT>( radius, angle, 3.0, 3.0 );
271 mx::improc::eigenImage<double> *noMask = nullptr;
272
273 const auto coordinates = mx::improc::annulusCoords<angleT>( radius, angle, 3.0, 3.0, 1.5, 2.5, 0.0, 360.0 );
274 const auto bufferedCoordinates =
275 mx::improc::annulusCoords<angleT>( radius, angle, 3.0, 3.0, 1.5, 2.5, 0.0, 360.0, noMask, 0.5 );
276
277 REQUIRE( coordinates.size() == 12 );
278 REQUIRE( bufferedCoordinates.size() == 24 );
279
280 const auto hasCoordinate = []( const std::vector<std::vector<int>> &region, int row, int column )
281 { return std::find( region.begin(), region.end(), std::vector<int>{ row, column } ) != region.end(); };
282
283 REQUIRE_FALSE( hasCoordinate( coordinates, 2, 3 ) );
284 REQUIRE( hasCoordinate( bufferedCoordinates, 2, 3 ) );
285 REQUIRE_FALSE( hasCoordinate( coordinates, 1, 1 ) );
286 REQUIRE( hasCoordinate( bufferedCoordinates, 1, 1 ) );
287 REQUIRE_FALSE( hasCoordinate( bufferedCoordinates, 0, 3 ) );
288
289 mx::improc::eigenImage<double> mask( 7, 7 );
290 mask.setOnes();
291 mask( 2, 3 ) = 0;
292 const auto maskedBufferedCoordinates =
293 mx::improc::annulusCoords<angleT>( radius, angle, 3.0, 3.0, 1.5, 2.5, 0.0, 360.0, &mask, 0.5 );
294 REQUIRE_FALSE( hasCoordinate( maskedBufferedCoordinates, 2, 3 ) );
295
296 const auto indices = mx::improc::annulusIndices<angleT>( radius, angle, 3.0, 3.0, 1.5, 2.5, 0.0, 360.0 );
297 const auto bufferedIndices =
298 mx::improc::annulusIndices<angleT>( radius, angle, 3.0, 3.0, 1.5, 2.5, 0.0, 360.0, noMask, 0.5 );
299
300 REQUIRE( indices.size() == coordinates.size() );
301 REQUIRE( bufferedIndices.size() == bufferedCoordinates.size() );
302 REQUIRE( std::find( indices.begin(), indices.end(), 23 ) == indices.end() );
303 REQUIRE( std::find( bufferedIndices.begin(), bufferedIndices.end(), 23 ) != bufferedIndices.end() );
304 REQUIRE( std::find( indices.begin(), indices.end(), 8 ) == indices.end() );
305 REQUIRE( std::find( bufferedIndices.begin(), bufferedIndices.end(), 8 ) != bufferedIndices.end() );
306
307 mx::improc::eigenImage<double> fractionalRadius( 10, 10 );
308 mx::improc::eigenImage<double> fractionalAngle;
309 mx::improc::radAngImage<angleT>( fractionalRadius, fractionalAngle, 4.5, 4.5 );
310
311 const auto fractionalCoordinates = mx::improc::annulusCoords<
312 angleT>( fractionalRadius, fractionalAngle, 4.5, 4.5, 3.0, 3.5, 0.0, 360.0, noMask, 0.5 );
313 REQUIRE( hasCoordinate( fractionalCoordinates, 1, 4 ) );
314 REQUIRE( hasCoordinate( fractionalCoordinates, 8, 4 ) );
315
316 const auto fractionalIndices = mx::improc::annulusIndices<
317 angleT>( fractionalRadius, fractionalAngle, 4.5, 4.5, 3.0, 3.5, 0.0, 360.0, noMask, 0.5 );
318 REQUIRE( std::find( fractionalIndices.begin(), fractionalIndices.end(), 41 ) != fractionalIndices.end() );
319 REQUIRE( std::find( fractionalIndices.begin(), fractionalIndices.end(), 48 ) != fractionalIndices.end() );
320
321 const auto emptyBufferedCoordinates =
322 mx::improc::annulusCoords<angleT>( radius, angle, 3.0, 3.0, 0.0, 0.5, 0.0, 360.0, noMask, -0.5 );
323 REQUIRE( emptyBufferedCoordinates.empty() );
324}
325
326/** \brief Verifies cutImageRegion and insertImageRegion preserve indexed pixel order and resize policy.
327 *
328 * \ingroup imageMasks_unit_tests
329 */
330TEST_CASE( "Image regions are cut and inserted by linear index", "[improc::imageMasks][imageRegion]" )
331{
332 mx::improc::eigenImage<double> image( 3, 4 );
333 for( Eigen::Index index = 0; index < image.size(); ++index )
334 {
335 image( index ) = static_cast<double>( 10 + index );
336 }
337 const std::vector<size_t> indices{ 0, 5, 11 };
338
339 mx::improc::eigenImage<double> cut;
340 mx::improc::cutImageRegion( cut, image, indices );
341 REQUIRE( cut.rows() == 3 );
342 REQUIRE( cut.cols() == 1 );
343 REQUIRE( cut( 0 ) == Approx( 10 ) );
344 REQUIRE( cut( 1 ) == Approx( 15 ) );
345 REQUIRE( cut( 2 ) == Approx( 21 ) );
346
347 mx::improc::eigenImage<double> preallocated( 3, 1 );
348 preallocated.setConstant( -1 );
349 mx::improc::cutImageRegion( preallocated, image, indices, false );
350 REQUIRE( preallocated.isApprox( cut ) );
351
352 mx::improc::eigenImage<double> inserted( image.size(), 1 );
353 inserted.setConstant( -1 );
354 auto insertedView = inserted.block( 0, 0, inserted.rows(), 1 );
355 mx::improc::insertImageRegion( insertedView, cut, indices );
356 REQUIRE( inserted( 0 ) == Approx( 10 ) );
357 REQUIRE( inserted( 5 ) == Approx( 15 ) );
358 REQUIRE( inserted( 11 ) == Approx( 21 ) );
359 REQUIRE( inserted( 1 ) == Approx( -1 ) );
360
361 mx::improc::cutImageRegion( cut, image, {} );
362 REQUIRE( cut.rows() == 0 );
363 REQUIRE( cut.cols() == 1 );
364 mx::improc::insertImageRegion( insertedView, cut, {} );
365 REQUIRE( inserted( 0 ) == Approx( 10 ) );
366}
eigenImage.hpp
Tools for using the eigen library for image processing.
mx::improc::eigenImage
Eigen::Array< scalarT, -1, -1 > eigenImage
Definition of the eigenImage type, which is an alias for Eigen::Array.
Definition eigenImage.hpp:44
mx::math::half_pi
constexpr T half_pi()
Get the value of pi/2.
Definition constants.hpp:152
TEST_CASE
TEST_CASE("Masking wedges in an image", "[improc::imageMasks::maskWedge]")
Definition imageMasks_test.cpp:29
mx::improc::cutImageRegion
void cutImageRegion(imageTout &imout, const imageTin &imin, const std::vector< size_t > &idx, bool resize=true)
Cut out a region of an image specified by an index-mask.
Definition imageMasks.hpp:972
mx::improc::insertImageRegion
void insertImageRegion(imageTout imout, const imageTin &imin, const std::vector< size_t > &idx)
Insert a region of an image specified by an index-mask.
Definition imageMasks.hpp:999
mx::improc::maskWedge
void maskWedge(arrayT &m, typename arrayT::Scalar xcen, typename arrayT::Scalar ycen, typename arrayT::Scalar angCen, typename arrayT::Scalar angHW, typename arrayT::Scalar val=0)
Mask a wedge in an image.
Definition imageMasks.hpp:681
mx::improc::annulusCoords
std::vector< std::vector< int > > annulusCoords(const eigenT1 &rIm, const eigenT2 &qIm, typename angleT::realT xcen, typename angleT::realT ycen, typename angleT::realT min_r, typename angleT::realT max_r, typename angleT::realT min_q, typename angleT::realT max_q, eigenT3 *mask=0, typename angleT::realT pixbuf=0)
Get the array coordinates of an annular region in an image.
Definition imageMasks.hpp:323
mx::improc::radAngImage
void radAngImage(eigenT1 &rIm, eigenT2 &qIm, typename angleT::realT xc, typename angleT::realT yc, typename angleT::realT rscale=1)
Fills in the cells of Eigen-like arrays with their radius amd angle relative to the center.
Definition imageMasks.hpp:152
mx::improc::annulusIndices
std::vector< size_t > annulusIndices(const eigenT1 &rIm, const eigenT2 &qIm, typename angleT::realT xcen, typename angleT::realT ycen, typename angleT::realT min_r, typename angleT::realT max_r, typename angleT::realT min_q, typename angleT::realT max_q, eigenT3 *mask=0, typename angleT::realT pixbuf=0)
Get the vector indices of an annular region in an image.
Definition imageMasks.hpp:368
imageMasks.hpp
Declares and defines functions to work with image masks.
mx::improc::rotateMask
void rotateMask(imageT &rotMask, imageT &mask, typename imageT::Scalar angle)
Rotate a binary mask.
Definition imageMasks.hpp:1017
mx::improc::cubicConvolTransform
Transformation by cubic convolution interpolation.
Definition imageTransforms.hpp:92