LSST Applications g180d380827+78227d2bc4,g2079a07aa2+86d27d4dc4,g2305ad1205+bdd7851fe3,g2bbee38e9b+c6a8a0fb72,g337abbeb29+c6a8a0fb72,g33d1c0ed96+c6a8a0fb72,g3a166c0a6a+c6a8a0fb72,g3d1719c13e+260d7c3927,g3ddfee87b4+723a6db5f3,g487adcacf7+29e55ea757,g50ff169b8f+96c6868917,g52b1c1532d+585e252eca,g591dd9f2cf+9443c4b912,g62aa8f1a4b+7e2ea9cd42,g858d7b2824+260d7c3927,g864b0138d7+8498d97249,g95921f966b+dffe86973d,g991b906543+260d7c3927,g99cad8db69+4809d78dd9,g9c22b2923f+e2510deafe,g9ddcbc5298+9a081db1e4,ga1e77700b3+03d07e1c1f,gb0e22166c9+60f28cb32d,gb23b769143+260d7c3927,gba4ed39666+c2a2e4ac27,gbb8dafda3b+e22341fd87,gbd998247f1+585e252eca,gc120e1dc64+713f94b854,gc28159a63d+c6a8a0fb72,gc3e9b769f7+385ea95214,gcf0d15dbbd+723a6db5f3,gdaeeff99f8+f9a426f77a,ge6526c86ff+fde82a80b9,ge79ae78c31+c6a8a0fb72,gee10cc3b42+585e252eca,w.2024.18
LSST Data Management Base Package
Loading...
Searching...
No Matches
Classes | Typedefs | Functions
detect_pybind11.cc File Reference
#include <pybind11/pybind11.h>
#include <pybind11/numpy.h>
#include <pybind11/stl.h>
#include <pybind11/eigen.h>
#include <math.h>
#include <algorithm>

Go to the source code of this file.

Classes

class  Peak
 A Peak in a Footprint This class is meant to keep track of both the position and flux at the location of a maximum in a Footprint. More...
 
class  Footprint
 

Typedefs

typedef Eigen::Array< int, 4, 1 > Bounds
 
typedef Eigen::Matrix< bool, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > MatrixB
 

Functions

template<typename M >
void get_connected_pixels (const int i, const int j, Eigen::Ref< const M > image, Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > unchecked, Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > footprint, Eigen::Ref< Bounds, 0, Eigen::Stride< 4, 1 > > bounds, const double thresh=0)
 
template<typename M >
MatrixB get_connected_multipeak (Eigen::Ref< const M > image, const std::vector< std::vector< int > > centers, const double thresh=0)
 Proximal operator to trim pixels not connected to one of the source centers.
 
bool sortBrightness (Peak a, Peak b)
 Sort two peaks, placing the brightest peak first.
 
template<typename M >
std::vector< Peakget_peaks (M &image, const double min_separation, const int y0, const int x0)
 
template<typename M >
void maskImage (Eigen::Ref< M, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > image, Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > footprint)
 
template<typename M , typename P >
std::vector< Footprintget_footprints (Eigen::Ref< const M > image, const double min_separation, const int min_area, const double thresh, const bool find_peaks=true)
 
 PYBIND11_MODULE (detect_pybind11, mod)
 

Typedef Documentation

◆ Bounds

typedef Eigen::Array<int, 4, 1> Bounds

Definition at line 11 of file detect_pybind11.cc.

◆ MatrixB

typedef Eigen::Matrix<bool, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> MatrixB

Definition at line 12 of file detect_pybind11.cc.

Function Documentation

◆ get_connected_multipeak()

template<typename M >
MatrixB get_connected_multipeak ( Eigen::Ref< const M > image,
const std::vector< std::vector< int > > centers,
const double thresh = 0 )

Proximal operator to trim pixels not connected to one of the source centers.

Definition at line 64 of file detect_pybind11.cc.

68 {
69 const int height = image.rows();
70 const int width = image.cols();
71 MatrixB unchecked = MatrixB::Ones(height, width);
72 MatrixB footprint = MatrixB::Zero(height, width);
73
74 for(auto center=begin(centers); center!=end(centers); ++center){
75 const int y = (*center)[0];
76 const int x = (*center)[1];
77 Bounds bounds; bounds << y, y, x, x;
78 get_connected_pixels(y, x, image, unchecked, footprint, bounds, thresh);
79 }
80
81 return footprint;
82}
int end
int y
Definition SpanSet.cc:48
T begin(T... args)
Eigen::Array< int, 4, 1 > Bounds
Eigen::Matrix< bool, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > MatrixB
void get_connected_pixels(const int i, const int j, Eigen::Ref< const M > image, Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > unchecked, Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > footprint, Eigen::Ref< Bounds, 0, Eigen::Stride< 4, 1 > > bounds, const double thresh=0)

◆ get_connected_pixels()

template<typename M >
void get_connected_pixels ( const int i,
const int j,
Eigen::Ref< const M > image,
Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > unchecked,
Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > footprint,
Eigen::Ref< Bounds, 0, Eigen::Stride< 4, 1 > > bounds,
const double thresh = 0 )

Definition at line 18 of file detect_pybind11.cc.

26 {
27 if(not unchecked(i,j)){
28 return;
29 }
30 unchecked(i,j) = false;
31
32 if(image(i,j) > thresh){
33 footprint(i,j) = true;
34
35 if(i < bounds[0]){
36 bounds[0] = i;
37 } else if(i > bounds[1]){
38 bounds[1] = i;
39 }
40 if(j < bounds[2]){
41 bounds[2] = j;
42 } else if(j > bounds[3]){
43 bounds[3] = j;
44 }
45
46 if(i > 0){
47 get_connected_pixels(i-1, j, image, unchecked, footprint, bounds, thresh);
48 }
49 if(i < image.rows()-1){
50 get_connected_pixels(i+1, j, image, unchecked, footprint, bounds, thresh);
51 }
52 if(j > 0){
53 get_connected_pixels(i, j-1, image, unchecked, footprint, bounds, thresh);
54 }
55 if(j < image.cols()-1){
56 get_connected_pixels(i, j+1, image, unchecked, footprint, bounds, thresh);
57 }
58 }
59}
afw::table::Key< afw::table::Array< ImagePixelT > > image

◆ get_footprints()

template<typename M , typename P >
std::vector< Footprint > get_footprints ( Eigen::Ref< const M > image,
const double min_separation,
const int min_area,
const double thresh,
const bool find_peaks = true )

Definition at line 246 of file detect_pybind11.cc.

252 {
253 const int height = image.rows();
254 const int width = image.cols();
255
256 std::vector<Footprint> footprints;
257 MatrixB unchecked = MatrixB::Ones(height, width);
258 MatrixB footprint = MatrixB::Zero(height, width);
259
260 for(int i=0; i<height; i++){
261 for(int j=0; j<width; j++){
262 Bounds bounds; bounds << i, i, j, j;
263 get_connected_pixels(i, j, image, unchecked, footprint, bounds, thresh);
264 int subHeight = bounds[1]-bounds[0]+1;
265 int subWidth = bounds[3]-bounds[2]+1;
266 if(subHeight * subWidth > min_area){
267 MatrixB subFootprint = footprint.block(bounds[0], bounds[2], subHeight, subWidth);
268 int area = subFootprint.count();
269 if(area >= min_area){
270 M patch = image.block(bounds[0], bounds[2], subHeight, subWidth);
271 maskImage<M>(patch, subFootprint);
272 std::vector<Peak> _peaks;
273 if(find_peaks){
274 _peaks = get_peaks(
275 patch,
276 min_separation,
277 bounds[0],
278 bounds[2]
279 );
280 }
281
282 footprints.push_back(Footprint(subFootprint, _peaks, bounds));
283 }
284 }
285 footprint.block(bounds[0], bounds[2], subHeight, subWidth) = MatrixB::Zero(subHeight, subWidth);
286 }
287 }
288 return footprints;
289}
std::vector< Peak > get_peaks(M &image, const double min_separation, const int y0, const int x0)

◆ get_peaks()

template<typename M >
std::vector< Peak > get_peaks ( M & image,
const double min_separation,
const int y0,
const int x0 )

Sort the peaks in the footprint so that the brightest are first

Definition at line 128 of file detect_pybind11.cc.

133 {
134 const int height = image.rows();
135 const int width = image.cols();
136
137 std::vector<Peak> peaks;
138
139 for(int i=0; i<height; i++){
140 for(int j=0; j<width; j++){
141 if(i > 0 && image(i, j) <= image(i-1, j)){
142 continue;
143 }
144 if(i < height-1 && image(i,j) <= image(i+1, j)){
145 continue;
146 }
147 if(j > 0 && image(i, j) <= image(i, j-1)){
148 continue;
149 }
150 if(j < width-1 && image(i,j) <= image(i, j+1)){
151 continue;
152 }
153
154 if(i > 0 && j > 0 && image(i, j) <= image(i-1, j-1)){
155 continue;
156 }
157 if(i < height-1 && j < width-1 && image(i,j) <= image(i+1, j+1)){
158 continue;
159 }
160 if(i < height-1 && j > 0 && image(i, j) <= image(i+1, j-1)){
161 continue;
162 }
163 if(i > 0 && j < width-1 && image(i,j) <= image(i-1, j+1)){
164 continue;
165 }
166
167 peaks.push_back(Peak(i+y0, j+x0, static_cast<double>(image(i, j))));
168 }
169 }
170
171 assert(peaks.size() > 0);
172
174 std::sort (peaks.begin(), peaks.end(), sortBrightness);
175
176 // Remove peaks within min_separation
177 double min_separation2 = min_separation * min_separation;
178 int i = 0;
179 while (i < peaks.size()-1){
180 int j = i+1;
181 Peak *p1 = &peaks[i];
182 while (j < peaks.size()){
183 Peak *p2 = &peaks[j];
184 double dy = p1->getY()-p2->getY();
185 double dx = p1->getX()-p2->getX();
186 double separation2 = dy*dy + dx*dx;
187 if(separation2 < min_separation2){
188 peaks.erase(peaks.begin()+j);
189 i--;
190 }
191 j++;
192 }
193 i++;
194 }
195
196 assert(peaks.size() > 0);
197
198 return peaks;
199}
A Peak in a Footprint This class is meant to keep track of both the position and flux at the location...
int getY()
bool sortBrightness(Peak a, Peak b)
Sort two peaks, placing the brightest peak first.
T end(T... args)
T erase(T... args)
T push_back(T... args)
T size(T... args)
T sort(T... args)

◆ maskImage()

template<typename M >
void maskImage ( Eigen::Ref< M, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > image,
Eigen::Ref< MatrixB, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > footprint )

Definition at line 228 of file detect_pybind11.cc.

231 {
232 const int height = image.rows();
233 const int width = image.cols();
234
235 for(int i=0; i<height; i++){
236 for(int j=0; j<width; j++){
237 if(!footprint(i,j)){
238 image(i,j) = 0;
239 }
240 }
241 }
242}

◆ PYBIND11_MODULE()

PYBIND11_MODULE ( detect_pybind11 ,
mod  )

Definition at line 293 of file detect_pybind11.cc.

293 {
294 typedef Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> MatrixF;
295 typedef Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> MatrixD;
296
297 mod.doc() = "Fast detection algorithms implemented in C++";
298
299 mod.def("get_connected_pixels", &get_connected_pixels<MatrixB>,
300 "Create a boolean mask for pixels that are connected");
301 mod.def("get_connected_pixels", &get_connected_pixels<MatrixF>,
302 "Create a boolean mask for pixels that are connected");
303 mod.def("get_connected_pixels", &get_connected_pixels<MatrixD>,
304 "Create a boolean mask for pixels that are connected");
305
306 mod.def("get_connected_multipeak", &get_connected_multipeak<MatrixB>,
307 "Trim pixels not conencted to a center from a list of centers");
308 mod.def("get_connected_multipeak", &get_connected_multipeak<MatrixF>,
309 "Trim pixels not conencted to a center from a list of centers");
310 mod.def("get_connected_multipeak", &get_connected_multipeak<MatrixD>,
311 "Trim pixels not conencted to a center from a list of centers");
312
313 mod.def("get_peaks", &get_peaks<MatrixF>,
314 "Get a list of peaks in a footprint created by get_connected_pixels");
315 mod.def("get_peaks", &get_peaks<MatrixD>,
316 "Get a list of peaks in a footprint created by get_connected_pixels");
317
318 mod.def("get_footprints", &get_footprints<MatrixF, float>,
319 "Create a list of all of the footprints in an image, with their peaks"
320 "image"_a, "min_separation"_a, "min_area"_a, "thresh"_a, "find_peaks"_a);
321 mod.def("get_footprints", &get_footprints<MatrixD, double>,
322 "Create a list of all of the footprints in an image, with their peaks"
323 "image"_a, "min_separation"_a, "min_area"_a, "thresh"_a, "find_peaks"_a);
324
325 py::class_<Footprint>(mod, "Footprint")
326 .def(py::init<MatrixB, std::vector<Peak>, Bounds>(),
327 "footprint"_a, "peaks"_a, "bounds"_a)
328 .def_property_readonly("data", &Footprint::getFootprint)
329 .def_readwrite("peaks", &Footprint::peaks)
330 .def_property_readonly("bounds", &Footprint::getBounds);
331
332 py::class_<Peak>(mod, "Peak")
333 .def(py::init<int, int, double>(),
334 "y"_a, "x"_a, "flux"_a)
335 .def_property_readonly("y", &Peak::getY)
336 .def_property_readonly("x", &Peak::getX)
337 .def_property_readonly("flux", &Peak::getFlux);
338}
std::vector< Peak > peaks
Bounds getBounds()
MatrixB getFootprint()
double getFlux()

◆ sortBrightness()

bool sortBrightness ( Peak a,
Peak b )

Sort two peaks, placing the brightest peak first.

Definition at line 117 of file detect_pybind11.cc.

117 {
118 return a.getFlux() > b.getFlux();
119}
table::Key< int > b