SourceXtractorPlusPlus  0.15
Please provide a description of the project.
GrowthCurveTask.cpp
Go to the documentation of this file.
1 
27 #include "SEUtils/Mat22.h"
28 
29 
30 namespace SourceXtractor {
31 
32 using SExtractor::Mat22;
33 
34 static const SeFloat GROWTH_NSIG = 6.;
35 static const size_t GROWTH_NSAMPLES = 64;
36 
37 static SeFloat getPixelValue(int x, int y, SeFloat centroid_x, SeFloat centroid_y,
39  const std::shared_ptr<ImageAccessor<SeFloat>>& variance_map, SeFloat variance_threshold,
40  bool use_symmetry) {
41  // Get the pixel value
42  DetectionImage::PixelType pixel_value = 0;
43  // Masked out
44  if (variance_map->getValue(x, y) > variance_threshold) {
45  if (use_symmetry) {
46  auto mirror_x = 2 * centroid_x - x + 0.49999;
47  auto mirror_y = 2 * centroid_y - y + 0.49999;
48  if (mirror_x >= 0 && mirror_y >= 0 && mirror_x < image->getWidth() && mirror_y < image->getHeight()) {
49  if (variance_map->getValue(mirror_x, mirror_y) < variance_threshold) {
50  // mirror pixel is OK: take the value
51  pixel_value = image->getValue(mirror_x, mirror_y);
52  }
53  }
54  }
55  }
56  // Not masked
57  else {
58  pixel_value = image->getValue(x, y);
59  }
60  return pixel_value;
61 }
62 
63 GrowthCurveTask::GrowthCurveTask(unsigned instance, bool use_symmetry)
64  : m_instance{instance}, m_use_symmetry{use_symmetry} {}
65 
67  const auto& measurement_frame_info = source.getProperty<MeasurementFrameInfo>(m_instance);
68  const auto& measurement_frame_images = source.getProperty<MeasurementFrameImages>(m_instance);
69 
70  auto variance_threshold = measurement_frame_info.getVarianceThreshold();
71 
72  const auto image = measurement_frame_images.getLockedImage(LayerSubtractedImage);
73  const auto variance_map = measurement_frame_images.getLockedImage(LayerVarianceMap);
74 
75  auto centroid_x = source.getProperty<MeasurementFramePixelCentroid>(m_instance).getCentroidX();
76  auto centroid_y = source.getProperty<MeasurementFramePixelCentroid>(m_instance).getCentroidY();
77  Mat22 jacobian{source.getProperty<JacobianSource>(m_instance).asTuple()};
78 
79  // These are from the detection frame
80  auto& shape = source.getProperty<ShapeParameters>();
81  auto& kron = source.getProperty<KronRadius>();
82 
83  // Radius for computing the growth curve and step size *on the detection frame*
84  double detection_rlim = std::max(GROWTH_NSIG * shape.getEllipseA(), kron.getKronRadius());
85 
86  // Now we need to compute the rlim for the measurement frame
87  // We take two vectors defined by the radius on the detection frame along the X and Y,
88  // transform them, and we use as a limit now the longest of the two after the transformation
89  Mat22 radius_22{detection_rlim, 0, 0, detection_rlim};
90  radius_22 = radius_22 * jacobian;
91  double r1 = radius_22[0] * radius_22[0] + radius_22[1] * radius_22[1];
92  double r2 = radius_22[2] * radius_22[2] + radius_22[3] * radius_22[3];
93  double rlim = std::sqrt(std::max(r1, r2));
94 
95  double step_size = rlim / GROWTH_NSAMPLES;
96 
97  // List of apertures
100  apertures.reserve(GROWTH_NSAMPLES);
101  for (size_t step = 1; step <= GROWTH_NSAMPLES; ++step) {
102  apertures.emplace_back(step_size * step);
103  }
104 
105  // Boundaries for the computation
106  // We know the last aperture is the widest, so we take the limits from it
107  auto min_coord = apertures.back().getMinPixel(centroid_x, centroid_y);
108  auto max_coord = apertures.back().getMaxPixel(centroid_x, centroid_y);
109 
110  // Compute fluxes for each ring
111  for (auto y = min_coord.m_y; y <= max_coord.m_y; ++y) {
112  for (auto x = min_coord.m_x; x <= max_coord.m_x; ++x) {
113  if (!image->isInside(x, y)) {
114  continue;
115  }
116 
117  auto pixel_value = getPixelValue(x, y, centroid_x, centroid_y, image,
118  variance_map, variance_threshold,
120 
121  // Assign the pixel value according to the affected area
122  auto dx = x - centroid_x;
123  auto dy = y - centroid_y;
124  double r = std::sqrt(dx * dx + dy * dy);
125 
126  // The pixel may be affected by multiple rings, so we look for those
127  // that overlap the start and the end of the pixels (which has size sqrt(2) on the diagonal)
128  size_t idx = 0;
129  if (r > 1.42 / 2.) {
130  idx = static_cast<size_t>((r - 1.42 / 2.) / step_size);
131  }
132  size_t outer_idx = static_cast<size_t>(std::ceil((r + 1.42 / 2.) / step_size));
133 
134  double inner = 0;
135  outer_idx = std::min(outer_idx, apertures.size() - 1);
136  for (; idx <= outer_idx; ++idx) {
137  auto& aperture = apertures[idx];
138  auto area = aperture.getArea(centroid_x, centroid_y, x, y);
139 
140  fluxes[idx] += area * pixel_value - inner;
141  inner = area * pixel_value;
142  }
143  }
144  }
145 
146  // Accumulate
147  for (size_t i = 1; i < fluxes.size(); ++i) {
148  fluxes[i] += fluxes[i - 1];
149  }
150 
151  // Set property
152  source.setIndexedProperty<GrowthCurve>(m_instance, std::move(fluxes), rlim);
153 }
154 
155 } // end of namespace SourceXtractor
MeasurementFrameInfo.h
Jacobian.h
SourceXtractor::ImageAccessor
Definition: ImageAccessor.h:41
std::shared_ptr
STL class.
std::move
T move(T... args)
SourceXtractor::Image< SeFloat >::PixelType
SeFloat PixelType
Definition: Image.h:47
MeasurementFramePixelCentroid.h
SourceXtractor::MeasurementFramePixelCentroid
Definition: MeasurementFramePixelCentroid.h:31
std::vector::reserve
T reserve(T... args)
SourceXtractor::SeFloat
SeFloat32 SeFloat
Definition: Types.h:32
std::vector
STL class.
std::vector::size
T size(T... args)
SourceXtractor::KronRadius
Kron radius.
Definition: KronRadius.h:36
SourceXtractor::GrowthCurveTask::computeProperties
void computeProperties(SourceInterface &source) const override
Computes one or more properties for the Source.
Definition: GrowthCurveTask.cpp:66
Mat22.h
std::vector::back
T back(T... args)
SourceXtractor::LayerVarianceMap
@ LayerVarianceMap
Definition: Frame.h:44
SourceXtractor::MeasurementFrameImages
Definition: MeasurementFrameImages.h:31
std::sqrt
T sqrt(T... args)
CircularAperture.h
SourceXtractor::GrowthCurveTask::m_instance
unsigned m_instance
Definition: GrowthCurveTask.h:34
SourceXtractor::getPixelValue
static SeFloat getPixelValue(int x, int y, SeFloat centroid_x, SeFloat centroid_y, const std::shared_ptr< ImageAccessor< SeFloat >> &image, const std::shared_ptr< ImageAccessor< SeFloat >> &variance_map, SeFloat variance_threshold, bool use_symmetry)
Definition: GrowthCurveTask.cpp:37
SourceXtractor
Definition: Aperture.h:30
dy
std::shared_ptr< EngineParameter > dy
Definition: MoffatModelFittingTask.cpp:93
SourceXtractor::GROWTH_NSIG
static const SeFloat GROWTH_NSIG
Definition: GrowthCurveTask.cpp:34
SourceXtractor::MeasurementFrameInfo
Definition: MeasurementFrameInfo.h:28
SourceXtractor::SourceInterface::setIndexedProperty
void setIndexedProperty(std::size_t index, Args... args)
Convenience template method to call setProperty() with a more user-friendly syntax.
Definition: SourceInterface.h:64
SourceXtractor::GROWTH_NSAMPLES
static const std::string GROWTH_NSAMPLES
Definition: GrowthCurveConfig.cpp:26
std::ceil
T ceil(T... args)
SourceXtractor::JacobianSource
Definition: Jacobian.h:51
std::min
T min(T... args)
SourceXtractor::GrowthCurveTask::GrowthCurveTask
GrowthCurveTask(unsigned instance, bool use_symmetry)
Definition: GrowthCurveTask.cpp:63
std::vector::emplace_back
T emplace_back(T... args)
GrowthCurve.h
SourceXtractor::GrowthCurve
Definition: GrowthCurve.h:30
x
std::shared_ptr< DependentParameter< std::shared_ptr< EngineParameter > > > x
Definition: MoffatModelFittingTask.cpp:94
GrowthCurveTask.h
SourceXtractor::SourceInterface::getProperty
const PropertyType & getProperty(unsigned int index=0) const
Convenience template method to call getProperty() with a more user-friendly syntax.
Definition: SourceInterface.h:57
SourceXtractor::LayerSubtractedImage
@ LayerSubtractedImage
Definition: Frame.h:38
ShapeParameters.h
SourceXtractor::SourceInterface
The SourceInterface is an abstract "source" that has properties attached to it.
Definition: SourceInterface.h:46
std::max
T max(T... args)
dx
std::shared_ptr< EngineParameter > dx
Definition: MoffatModelFittingTask.cpp:93
SourceXtractor::GrowthCurveTask::m_use_symmetry
bool m_use_symmetry
Definition: GrowthCurveTask.h:35
y
std::shared_ptr< DependentParameter< std::shared_ptr< EngineParameter > > > y
Definition: MoffatModelFittingTask.cpp:94
SourceXtractor::ShapeParameters
Definition: ShapeParameters.h:32
MeasurementFrameImages.h
SExtractor::Mat22
Definition: Mat22.h:19
KronRadius.h