LSST Applications  21.0.0+04719a4bac,21.0.0-1-ga51b5d4+f5e6047307,21.0.0-11-g2b59f77+a9c1acf22d,21.0.0-11-ga42c5b2+86977b0b17,21.0.0-12-gf4ce030+76814010d2,21.0.0-13-g1721dae+760e7a6536,21.0.0-13-g3a573fe+768d78a30a,21.0.0-15-g5a7caf0+f21cbc5713,21.0.0-16-g0fb55c1+b60e2d390c,21.0.0-19-g4cded4ca+71a93a33c0,21.0.0-2-g103fe59+bb20972958,21.0.0-2-g45278ab+04719a4bac,21.0.0-2-g5242d73+3ad5d60fb1,21.0.0-2-g7f82c8f+8babb168e8,21.0.0-2-g8f08a60+06509c8b61,21.0.0-2-g8faa9b5+616205b9df,21.0.0-2-ga326454+8babb168e8,21.0.0-2-gde069b7+5e4aea9c2f,21.0.0-2-gecfae73+1d3a86e577,21.0.0-2-gfc62afb+3ad5d60fb1,21.0.0-25-g1d57be3cd+e73869a214,21.0.0-3-g357aad2+ed88757d29,21.0.0-3-g4a4ce7f+3ad5d60fb1,21.0.0-3-g4be5c26+3ad5d60fb1,21.0.0-3-g65f322c+e0b24896a3,21.0.0-3-g7d9da8d+616205b9df,21.0.0-3-ge02ed75+a9c1acf22d,21.0.0-4-g591bb35+a9c1acf22d,21.0.0-4-g65b4814+b60e2d390c,21.0.0-4-gccdca77+0de219a2bc,21.0.0-4-ge8a399c+6c55c39e83,21.0.0-5-gd00fb1e+05fce91b99,21.0.0-6-gc675373+3ad5d60fb1,21.0.0-64-g1122c245+4fb2b8f86e,21.0.0-7-g04766d7+cd19d05db2,21.0.0-7-gdf92d54+04719a4bac,21.0.0-8-g5674e7b+d1bd76f71f,master-gac4afde19b+a9c1acf22d,w.2021.13
LSST Data Management Base Package
FootprintMerge.cc
Go to the documentation of this file.
1 /*
2  * LSST Data Management System
3  * Copyright 2008-2014 LSST Corporation.
4  *
5  * This product includes software developed by the
6  * LSST Project (http://www.lsst.org/).
7  *
8  * This program is free software: you can redistribute it and/or modify
9  * it under the terms of the GNU General Public License as published by
10  * the Free Software Foundation, either version 3 of the License, or
11  * (at your option) any later version.
12  *
13  * This program is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the LSST License Statement and
19  * the GNU General Public License along with this program. If not,
20  * see <http://www.lsstcorp.org/LegalNotices/>.
21  */
22 #include <cstdint>
23 
24 #include "boost/bind.hpp"
25 
29 
30 namespace lsst {
31 namespace afw {
32 namespace detection {
33 
35 public:
36  typedef FootprintMergeList::KeyTuple KeyTuple;
38 
42  afw::table::SchemaMapper const &peakSchemaMapper, KeyTuple const &keys)
43  : _footprints(1, footprint), _source(sourceTable->makeRecord()) {
44  std::shared_ptr<Footprint> newFootprint = std::make_shared<Footprint>(*footprint);
45 
46  _source->set(keys.footprint, true);
47  // Replace all the Peaks in the merged Footprint with new ones that include the origin flags
48  newFootprint->getPeaks() = PeakCatalog(peakTable);
49  for (PeakCatalog::iterator iter = footprint->getPeaks().begin(); iter != footprint->getPeaks().end();
50  ++iter) {
51  std::shared_ptr<PeakRecord> newPeak = peakTable->copyRecord(*iter, peakSchemaMapper);
52  newPeak->set(keys.peak, true);
53  newFootprint->getPeaks().push_back(newPeak);
54  }
55  _source->setFootprint(newFootprint);
56  }
57 
58  ~FootprintMerge() = default;
59  FootprintMerge(FootprintMerge const &) = default;
63 
64  /*
65  * Does this Footprint overlap the merged Footprint.
66  *
67  * The current implementation just builds an image from the two Footprints and
68  * detects the number of peaks. This is not very efficient and will be changed
69  * within the Footprint class in the future.
70  */
71  bool overlaps(Footprint const &rhs) const {
72  return getMergedFootprint()->getSpans()->overlaps(*(rhs.getSpans()));
73  }
74 
75  /*
76  * Add this Footprint to the merge.
77  *
78  * If minNewPeakDist >= 0, it will add all peaks from the footprint to the merged Footprint
79  * that are greater than minNewPeakDist away from the closest peak in the existing list.
80  * If minNewPeakDist < 0, no peaks will be added from the footprint.
81  *
82  * If maxSamePeakDist >= 0, it will find the closest peak in the existing list for every peak
83  * in the footprint. If the closest peak is less than maxSamePeakDist, the peak will not
84  * be added and the closest peak will be flagged as detected by the filter defined in keys.
85  *
86  * If the footprint does not overlap it will do nothing.
87  */
88  void add(std::shared_ptr<Footprint> footprint, afw::table::SchemaMapper const &peakSchemaMapper,
89  KeyTuple const &keys, float minNewPeakDist = -1., float maxSamePeakDist = -1.) {
90  if (addSpans(footprint)) {
91  _footprints.push_back(footprint);
92  _source->set(keys.footprint, true);
93  _addPeaks(footprint->getPeaks(), &peakSchemaMapper, &keys, minNewPeakDist, maxSamePeakDist, NULL);
94  }
95  }
96 
97  /*
98  * Merge an already-merged clump of Footprints into this
99  *
100  * If minNewPeakDist >= 0, it will add all peaks from the footprint to the merged Footprint
101  * that are greater than minNewPeakDist away from the closest peak in the existing list.
102  * If minNewPeakDist < 0, no peaks will be added from the footprint.
103  *
104  * If maxSamePeakDist >= 0, it will find the closest peak in the existing list for every peak
105  * in the footprint. If the closest peak is less than maxSamePeakDist, the peak will not
106  * be added to the list and the flags from the closest peak will be set to the OR of the two.
107  *
108  * If the FootprintMerge does not overlap it will do nothing.
109  */
110  void add(FootprintMerge const &other, FilterMap const &keys, float minNewPeakDist = -1.,
111  float maxSamePeakDist = -1.) {
112  if (addSpans(other.getMergedFootprint())) {
113  _footprints.insert(_footprints.end(), other._footprints.begin(), other._footprints.end());
114  // Set source flags to the OR of the flags of the two inputs
115  for (FilterMap::const_iterator i = keys.begin(); i != keys.end(); ++i) {
116  afw::table::Key<afw::table::Flag> const &flagKey = i->second.footprint;
117  _source->set(flagKey, _source->get(flagKey) || other._source->get(flagKey));
118  }
119  _addPeaks(other.getMergedFootprint()->getPeaks(), NULL, NULL, minNewPeakDist, maxSamePeakDist,
120  &keys);
121  }
122  }
123 
124  // Get the bounding box of the merge
125  lsst::geom::Box2I getBBox() const { return getMergedFootprint()->getBBox(); }
126 
127  std::shared_ptr<Footprint> getMergedFootprint() const { return _source->getFootprint(); }
128 
130 
131  // Implementation helper for add() methods; returns true if the Footprint actually overlapped
132  // and was merged, and false otherwise.
134  if (!getMergedFootprint()->getSpans()->overlaps(*(footprint->getSpans()))) return false;
135  getMergedFootprint()->setSpans(getMergedFootprint()->getSpans()->union_(*(footprint->getSpans())));
136  return true;
137  }
138 
139 private:
140  /*
141  * Add new peaks to the list of peaks of the merged footprint.
142  * This function handles two different cases:
143  * - The peaks come from a single footprint. In this case, the peakSchemaMapper
144  * and keys should be defined so that it can create a new peak, copy the appropriate
145  * data, and set the peak flag defined in keys.
146  * - The peaks come from another FootprintMerge. In this case, filterMap should
147  * be defined so that the information from the other peaks can be propagated.
148  */
149  void _addPeaks(PeakCatalog const &otherPeaks, afw::table::SchemaMapper const *peakSchemaMapper,
150  KeyTuple const *keys, float minNewPeakDist, float maxSamePeakDist,
151  FilterMap const *filterMap) {
152  if (minNewPeakDist < 0 && maxSamePeakDist < 0) return;
153 
154  assert(peakSchemaMapper || filterMap);
155 
156  PeakCatalog &currentPeaks = getMergedFootprint()->getPeaks();
157  std::shared_ptr<PeakRecord> nearestPeak;
158  // Create new list of peaks
159  PeakCatalog newPeaks(currentPeaks.getTable());
160  float minNewPeakDist2 = minNewPeakDist * minNewPeakDist;
161  float maxSamePeakDist2 = maxSamePeakDist * maxSamePeakDist;
162  for (PeakCatalog::const_iterator otherIter = otherPeaks.begin(); otherIter != otherPeaks.end();
163  ++otherIter) {
164  float minDist2 = std::numeric_limits<float>::infinity();
165 
166  for (PeakCatalog::const_iterator currentIter = currentPeaks.begin();
167  currentIter != currentPeaks.end(); ++currentIter) {
168  float dist2 = otherIter->getI().distanceSquared(currentIter->getI());
169 
170  if (dist2 < minDist2) {
171  minDist2 = dist2;
172  nearestPeak = currentIter;
173  }
174  }
175 
176  if (minDist2 < maxSamePeakDist2 && nearestPeak && maxSamePeakDist > 0) {
177  if (peakSchemaMapper) {
178  nearestPeak->set(keys->peak, true);
179  } else {
180  for (FilterMap::const_iterator i = filterMap->begin(); i != filterMap->end(); ++i) {
181  afw::table::Key<afw::table::Flag> const &flagKey = i->second.peak;
182  nearestPeak->set(flagKey, nearestPeak->get(flagKey) || otherIter->get(flagKey));
183  }
184  }
185  } else if (minDist2 > minNewPeakDist2 && !(minNewPeakDist < 0)) {
186  if (peakSchemaMapper) {
187  std::shared_ptr<PeakRecord> newPeak = newPeaks.addNew();
188  newPeak->assign(*otherIter, *peakSchemaMapper);
189  newPeak->set(keys->peak, true);
190  } else {
191  newPeaks.push_back(otherIter);
192  }
193  }
194  }
195 
196  getMergedFootprint()->getPeaks().insert(getMergedFootprint()->getPeaks().end(), newPeaks.begin(),
197  newPeaks.end(),
198  true // deep-copy
199  );
200  }
201 
204 };
205 
207  std::vector<std::string> const &filterList,
208  afw::table::Schema const &initialPeakSchema)
209  : _peakSchemaMapper(initialPeakSchema) {
210  _initialize(sourceSchema, filterList);
211 }
212 
214  std::vector<std::string> const &filterList)
215  : _peakSchemaMapper(PeakTable::makeMinimalSchema()) {
216  _initialize(sourceSchema, filterList);
217 }
218 
224 
225 void FootprintMergeList::_initialize(afw::table::Schema &sourceSchema,
226  std::vector<std::string> const &filterList) {
227  _peakSchemaMapper.addMinimalSchema(_peakSchemaMapper.getInputSchema(), true);
228  // Add Flags for the filters
229  for (std::vector<std::string>::const_iterator iter = filterList.begin(); iter != filterList.end();
230  ++iter) {
231  KeyTuple &keys = _filterMap[*iter];
232  keys.footprint = sourceSchema.addField<afw::table::Flag>(
233  "merge_footprint_" + *iter,
234  "Detection footprint overlapped with a detection from filter " + *iter);
235  keys.peak = _peakSchemaMapper.editOutputSchema().addField<afw::table::Flag>(
236  "merge_peak_" + *iter, "Peak detected in filter " + *iter);
237  }
238  _peakTable = PeakTable::make(_peakSchemaMapper.getOutputSchema());
239 }
240 
242  afw::table::SourceCatalog const &inputCat, std::string const &filter,
243  float minNewPeakDist, bool doMerge, float maxSamePeakDist) {
244  FilterMap::const_iterator keyIter = _filterMap.find(filter);
245  if (keyIter == _filterMap.end()) {
247  (boost::format("Filter %s not in original list") % filter).str());
248  }
249 
250  // If list is empty or merging not requested, don't check for any matches, just add all the objects
251  bool checkForMatches = !_mergeList.empty() && doMerge;
252 
253  for (afw::table::SourceCatalog::const_iterator srcIter = inputCat.begin(); srcIter != inputCat.end();
254  ++srcIter) {
255  // Only consider unblended objects
256  if (srcIter->getParent() != 0) continue;
257 
258  std::shared_ptr<Footprint> foot = srcIter->getFootprint();
259 
260  // Empty pointer to account for the first match in the catalog. If there is more than one
261  // match, subsequent matches will be merged with this one
263 
264  if (checkForMatches) {
265  FootprintMergeVec::iterator iter = _mergeList.begin();
266  while (iter != _mergeList.end()) {
267  // Grow by one pixel to allow for touching
268  lsst::geom::Box2I box((**iter).getBBox());
269  box.grow(lsst::geom::Extent2I(1, 1));
270  if (box.overlaps(foot->getBBox()) && (**iter).overlaps(*foot)) {
271  if (!first) {
272  first = *iter;
273  // Spatially extend existing FootprintMerge in order to connect subsequent,
274  // now-overlapping FootprintMerges. If a subsequent FootprintMerge overlaps with
275  // the new footprint, it's now guaranteed to overlap with this first FootprintMerge.
276  // Hold off adding foot's lower-priority footprints and peaks until the
277  // higher-priority existing peaks are merged into this first FootprintMerge.
278  first->addSpans(foot);
279  } else {
280  // Add existing merged Footprint to first
281  first->add(**iter, _filterMap, minNewPeakDist, maxSamePeakDist);
282  iter = _mergeList.erase(iter);
283  continue;
284  }
285  }
286  ++iter;
287  } // while mergeList
288  } // if checkForMatches
289 
290  if (first) {
291  // Now merge footprint including peaks into the newly-connected, higher-priority FootprintMerge
292  first->add(foot, _peakSchemaMapper, keyIter->second, minNewPeakDist, maxSamePeakDist);
293  } else {
294  // Footprint did not overlap with any existing FootprintMerges. Add to MergeList
295  _mergeList.push_back(std::make_shared<FootprintMerge>(foot, sourceTable, _peakTable,
296  _peakSchemaMapper, keyIter->second));
297  }
298  }
299 }
300 
302  // Now set the merged footprint as the footprint of the SourceRecord
303  for (FootprintMergeVec::iterator iter = _mergeList.begin(); iter != _mergeList.end(); ++iter) {
304  outputCat.push_back((**iter).getSource());
305  }
306 }
307 } // namespace detection
308 } // namespace afw
309 } // namespace lsst
int end
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition: Exception.h:48
ItemVariant const * other
Definition: Schema.cc:56
T begin(T... args)
Class to describe the properties of a detected object from an image.
Definition: Footprint.h:63
std::shared_ptr< geom::SpanSet > getSpans() const
Return a shared pointer to the SpanSet.
Definition: Footprint.h:115
bool overlaps(Footprint const &rhs) const
bool addSpans(std::shared_ptr< Footprint > footprint)
std::shared_ptr< afw::table::SourceRecord > getSource() const
FootprintMergeList::KeyTuple KeyTuple
std::shared_ptr< Footprint > getMergedFootprint() const
FootprintMerge(FootprintMerge &&)=default
FootprintMerge(std::shared_ptr< Footprint > footprint, std::shared_ptr< afw::table::SourceTable > sourceTable, std::shared_ptr< PeakTable > peakTable, afw::table::SchemaMapper const &peakSchemaMapper, KeyTuple const &keys)
FootprintMerge & operator=(FootprintMerge &&)=default
void add(std::shared_ptr< Footprint > footprint, afw::table::SchemaMapper const &peakSchemaMapper, KeyTuple const &keys, float minNewPeakDist=-1., float maxSamePeakDist=-1.)
FootprintMergeList::FilterMap FilterMap
lsst::geom::Box2I getBBox() const
FootprintMerge(FootprintMerge const &)=default
FootprintMerge & operator=(FootprintMerge const &)=default
void add(FootprintMerge const &other, FilterMap const &keys, float minNewPeakDist=-1., float maxSamePeakDist=-1.)
FootprintMergeList(afw::table::Schema &sourceSchema, std::vector< std::string > const &filterList, afw::table::Schema const &initialPeakSchema)
Initialize the merge with a custom initial peak schema.
void getFinalSources(afw::table::SourceCatalog &outputCat)
Get SourceCatalog with entries that contain the final Footprint and SourceRecord for each entry.
FootprintMergeList & operator=(FootprintMergeList const &)
void addCatalog(std::shared_ptr< afw::table::SourceTable > sourceTable, afw::table::SourceCatalog const &inputCat, std::string const &filter, float minNewPeakDist=-1., bool doMerge=true, float maxSamePeakDist=-1.)
Add objects from a SourceCatalog in the specified filter.
Table class for Peaks in Footprints.
Definition: Peak.h:102
static std::shared_ptr< PeakTable > make(afw::table::Schema const &schema, bool forceNew=false)
Obtain a table that can be used to create records with given schema.
Definition: Peak.cc:98
Iterator class for CatalogT.
Definition: Catalog.h:39
iterator begin()
Iterator access.
Definition: Catalog.h:396
std::shared_ptr< Table > getTable() const
Return the table associated with the catalog.
Definition: Catalog.h:114
Defines the fields and offsets for a table.
Definition: Schema.h:50
Key< T > addField(Field< T > const &field, bool doReplace=false)
Add a new field to the Schema, and return the associated Key.
Definition: Schema.cc:668
A mapping between the keys of two Schemas, used to copy data between them.
Definition: SchemaMapper.h:21
Schema const getOutputSchema() const
Return the output schema (copy-on-write).
Definition: SchemaMapper.h:27
Schema & editOutputSchema()
Return a reference to the output schema that allows it to be modified in place.
Definition: SchemaMapper.h:30
void addMinimalSchema(Schema const &minimal, bool doMap=true)
Add the given minimal schema to the output schema.
Schema const getInputSchema() const
Return the input schema (copy-on-write).
Definition: SchemaMapper.h:24
An integer coordinate rectangle.
Definition: Box.h:55
bool overlaps(Box2I const &other) const noexcept
Return true if any points in other are also in this.
Definition: Box.cc:122
void grow(int buffer)
Increase the size of the box by the given buffer amount in all directions.
Definition: Box.h:249
Reports errors in the logical structure of the program.
Definition: Runtime.h:46
T empty(T... args)
T end(T... args)
T erase(T... args)
T find(T... args)
T get(T... args)
T infinity(T... args)
afw::table::CatalogT< PeakRecord > PeakCatalog
Definition: Peak.h:244
def format(config, name=None, writeSourceLine=True, prefix="", verbose=False)
Definition: history.py:174
A base class for image defects.
T push_back(T... args)