LSST Applications 27.0.0,g0265f82a02+469cd937ee,g02d81e74bb+21ad69e7e1,g1470d8bcf6+cbe83ee85a,g2079a07aa2+e67c6346a6,g212a7c68fe+04a9158687,g2305ad1205+94392ce272,g295015adf3+81dd352a9d,g2bbee38e9b+469cd937ee,g337abbeb29+469cd937ee,g3939d97d7f+72a9f7b576,g487adcacf7+71499e7cba,g50ff169b8f+5929b3527e,g52b1c1532d+a6fc98d2e7,g591dd9f2cf+df404f777f,g5a732f18d5+be83d3ecdb,g64a986408d+21ad69e7e1,g858d7b2824+21ad69e7e1,g8a8a8dda67+a6fc98d2e7,g99cad8db69+f62e5b0af5,g9ddcbc5298+d4bad12328,ga1e77700b3+9c366c4306,ga8c6da7877+71e4819109,gb0e22166c9+25ba2f69a1,gb6a65358fc+469cd937ee,gbb8dafda3b+69d3c0e320,gc07e1c2157+a98bf949bb,gc120e1dc64+615ec43309,gc28159a63d+469cd937ee,gcf0d15dbbd+72a9f7b576,gdaeeff99f8+a38ce5ea23,ge6526c86ff+3a7c1ac5f1,ge79ae78c31+469cd937ee,gee10cc3b42+a6fc98d2e7,gf1cff7945b+21ad69e7e1,gfbcc870c63+9a11dc8c8f
LSST Data Management Base Package
Loading...
Searching...
No Matches
_saturated.cc
Go to the documentation of this file.
1/*
2 * Handle saturated pixels when making colour images
3 */
4#include <cmath>
5#include "boost/format.hpp"
8#include "Rgb.h"
9
10namespace lsst {
11namespace afw {
12namespace display {
13
14namespace {
15template <typename ImageT>
16class SetPixels {
17public:
18 explicit SetPixels() {}
19
20 void setValue(float value) { _value = value; }
21
22 void operator()(lsst::geom::Point2I const& point, typename ImageT::Pixel& arrayInput) {
23 arrayInput = _value;
24 }
25
26private:
27 float _value{0};
28};
29} // namespace
30
31template <typename ImageT>
32void replaceSaturatedPixels(ImageT& rim, // R image (e.g. i)
33 ImageT& gim, // G image (e.g. r)
34 ImageT& bim, // B image (e.g. g)
35 int borderWidth, // width of border used to estimate colour of saturated regions
36 float saturatedPixelValue // the brightness of a saturated pixel, once fixed
37) {
38 int const width = rim.getWidth(), height = rim.getHeight();
39 int const x0 = rim.getX0(), y0 = rim.getY0();
40
41 if (width != gim.getWidth() || height != gim.getHeight() || x0 != gim.getX0() || y0 != gim.getY0()) {
42 throw LSST_EXCEPT(
44 str(boost::format("R image has different size/origin from G image "
45 "(%dx%d+%d+%d v. %dx%d+%d+%d") %
46 width % height % x0 % y0 % gim.getWidth() % gim.getHeight() % gim.getX0() % gim.getY0()));
47 }
48 if (width != bim.getWidth() || height != bim.getHeight() || x0 != bim.getX0() || y0 != bim.getY0()) {
49 throw LSST_EXCEPT(
51 str(boost::format("R image has different size/origin from B image "
52 "(%dx%d+%d+%d v. %dx%d+%d+%d") %
53 width % height % x0 % y0 % bim.getWidth() % bim.getHeight() % bim.getX0() % bim.getY0()));
54 }
55
56 bool const useMaxPixel = !std::isfinite(saturatedPixelValue);
57
58 SetPixels<typename ImageT::Image> setR, setG, setB; // functors used to set pixel values
59
60 // Find all the saturated pixels in any of the three image
61 int const npixMin = 1; // minimum number of pixels in an object
62 afw::image::MaskPixel const SAT = rim.getMask()->getPlaneBitMask("SAT");
64
65 detection::FootprintSet sat(*rim.getMask(), satThresh, npixMin);
66 sat.merge(detection::FootprintSet(*gim.getMask(), satThresh, npixMin));
67 sat.merge(detection::FootprintSet(*bim.getMask(), satThresh, npixMin));
68 // go through the list of saturated regions, determining the mean colour of the surrounding pixels
69 using FootprintList = detection::FootprintSet::FootprintList;
71 for (const auto& foot : *feet) {
72 auto const bigFoot = std::make_shared<detection::Footprint>(foot->getSpans()->dilated(borderWidth),
73 foot->getRegion());
74
75 double sumR = 0, sumG = 0, sumB = 0; // sum of all non-saturated adjoining pixels
76 double maxR = 0, maxG = 0, maxB = 0; // maximum of non-saturated adjoining pixels
77
78 for (auto span : *bigFoot->getSpans()) {
79 int const y = span.getY() - y0;
80 if (y < 0 || y >= height) {
81 continue;
82 }
83 int sx0 = span.getX0() - x0;
84 if (sx0 < 0) {
85 sx0 = 0;
86 }
87 int sx1 = span.getX1() - x0;
88 if (sx1 >= width) {
89 sx1 = width - 1;
90 }
91
92 for (typename ImageT::iterator rptr = rim.at(sx0, y), rend = rim.at(sx1 + 1, y),
93 gptr = gim.at(sx0, y), bptr = bim.at(sx0, y);
94 rptr != rend; ++rptr, ++gptr, ++bptr) {
95 if (!((rptr.mask() | gptr.mask() | bptr.mask()) & SAT)) {
96 float val = rptr.image();
97 sumR += val;
98 if (val > maxR) {
99 maxR = val;
100 }
101
102 val = gptr.image();
103 sumG += val;
104 if (val > maxG) {
105 maxG = val;
106 }
107
108 val = bptr.image();
109 sumB += val;
110 if (val > maxB) {
111 maxB = val;
112 }
113 }
114 }
115 }
116 // OK, we have the mean fluxes for the pixels surrounding this set of saturated pixels
117 // so we can figure out the proper values to use for the saturated ones
118 float R = 0, G = 0, B = 0; // mean intensities
119 if (sumR + sumB + sumG > 0) {
120 if (sumR > sumG) {
121 if (sumR > sumB) {
122 R = useMaxPixel ? maxR : saturatedPixelValue;
123
124 G = (R * sumG) / sumR;
125 B = (R * sumB) / sumR;
126 } else {
127 B = useMaxPixel ? maxB : saturatedPixelValue;
128 R = (B * sumR) / sumB;
129 G = (B * sumG) / sumB;
130 }
131 } else {
132 if (sumG > sumB) {
133 G = useMaxPixel ? maxG : saturatedPixelValue;
134 R = (G * sumR) / sumG;
135 B = (G * sumB) / sumG;
136 } else {
137 B = useMaxPixel ? maxB : saturatedPixelValue;
138 R = (B * sumR) / sumB;
139 G = (B * sumG) / sumB;
140 }
141 }
142 }
143 // Now that we know R, G, and B we can fix the values
144 setR.setValue(R);
145 foot->getSpans()->applyFunctor(setR, *rim.getImage());
146 setG.setValue(G);
147 foot->getSpans()->applyFunctor(setG, *gim.getImage());
148 setB.setValue(B);
149 foot->getSpans()->applyFunctor(setB, *bim.getImage());
150 }
151}
152
154 image::MaskedImage<float>& bim, int borderWidth,
155 float saturatedPixelValue);
156} // namespace display
157} // namespace afw
158} // namespace lsst
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition Exception.h:48
int y
Definition SpanSet.cc:48
A set of Footprints, associated with a MaskedImage.
std::shared_ptr< FootprintList > getFootprints()
: Return the Footprints of detected objects
void merge(FootprintSet const &rhs, int tGrow=0, int rGrow=0, bool isotropic=true)
Merge a FootprintSet into *this.
std::vector< std::shared_ptr< Footprint > > FootprintList
The FootprintSet's set of Footprints.
A Threshold is used to pass a threshold value to detection algorithms.
Definition Threshold.h:43
@ BITMASK
Use (pixels & (given mask))
Definition Threshold.h:48
A class to manipulate images, masks, and variance as a single object.
Definition MaskedImage.h:74
Reports invalid arguments.
Definition Runtime.h:66
T isfinite(T... args)
void replaceSaturatedPixels(ImageT &rim, ImageT &gim, ImageT &bim, int borderWidth, float saturatedPixelValue)
Definition _saturated.cc:32
ImageT val
Definition CR.cc:146