LSST Applications  21.0.0-147-g0e635eb1+1acddb5be5,22.0.0+052faf71bd,22.0.0+1ea9a8b2b2,22.0.0+6312710a6c,22.0.0+729191ecac,22.0.0+7589c3a021,22.0.0+9f079a9461,22.0.1-1-g7d6de66+b8044ec9de,22.0.1-1-g87000a6+536b1ee016,22.0.1-1-g8e32f31+6312710a6c,22.0.1-10-gd060f87+016f7cdc03,22.0.1-12-g9c3108e+df145f6f68,22.0.1-16-g314fa6d+c825727ab8,22.0.1-19-g93a5c75+d23f2fb6d8,22.0.1-19-gb93eaa13+aab3ef7709,22.0.1-2-g8ef0a89+b8044ec9de,22.0.1-2-g92698f7+9f079a9461,22.0.1-2-ga9b0f51+052faf71bd,22.0.1-2-gac51dbf+052faf71bd,22.0.1-2-gb66926d+6312710a6c,22.0.1-2-gcb770ba+09e3807989,22.0.1-20-g32debb5+b8044ec9de,22.0.1-23-gc2439a9a+fb0756638e,22.0.1-3-g496fd5d+09117f784f,22.0.1-3-g59f966b+1e6ba2c031,22.0.1-3-g849a1b8+f8b568069f,22.0.1-3-gaaec9c0+c5c846a8b1,22.0.1-32-g5ddfab5d3+60ce4897b0,22.0.1-4-g037fbe1+64e601228d,22.0.1-4-g8623105+b8044ec9de,22.0.1-5-g096abc9+d18c45d440,22.0.1-5-g15c806e+57f5c03693,22.0.1-7-gba73697+57f5c03693,master-g6e05de7fdc+c1283a92b8,master-g72cdda8301+729191ecac,w.2021.39
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 
26 
27 namespace lsst {
28 namespace afw {
29 namespace detection {
30 
32 public:
33  using KeyTuple = FootprintMergeList::KeyTuple;
35 
39  afw::table::SchemaMapper const &peakSchemaMapper, KeyTuple const &keys)
40  : _footprints(1, footprint), _source(sourceTable->makeRecord()) {
41  std::shared_ptr<Footprint> newFootprint = std::make_shared<Footprint>(*footprint);
42 
43  _source->set(keys.footprint, true);
44  // Replace all the Peaks in the merged Footprint with new ones that include the origin flags
45  newFootprint->getPeaks() = PeakCatalog(peakTable);
46  for (PeakCatalog::iterator iter = footprint->getPeaks().begin(); iter != footprint->getPeaks().end();
47  ++iter) {
48  std::shared_ptr<PeakRecord> newPeak = peakTable->copyRecord(*iter, peakSchemaMapper);
49  newPeak->set(keys.peak, true);
50  newFootprint->getPeaks().push_back(newPeak);
51  }
52  _source->setFootprint(newFootprint);
53  }
54 
55  ~FootprintMerge() = default;
56  FootprintMerge(FootprintMerge const &) = default;
60 
61  /*
62  * Does this Footprint overlap the merged Footprint.
63  *
64  * The current implementation just builds an image from the two Footprints and
65  * detects the number of peaks. This is not very efficient and will be changed
66  * within the Footprint class in the future.
67  */
68  bool overlaps(Footprint const &rhs) const {
69  return getMergedFootprint()->getSpans()->overlaps(*(rhs.getSpans()));
70  }
71 
72  /*
73  * Add this Footprint to the merge.
74  *
75  * If minNewPeakDist >= 0, it will add all peaks from the footprint to the merged Footprint
76  * that are greater than minNewPeakDist away from the closest peak in the existing list.
77  * If minNewPeakDist < 0, no peaks will be added from the footprint.
78  *
79  * If maxSamePeakDist >= 0, it will find the closest peak in the existing list for every peak
80  * in the footprint. If the closest peak is less than maxSamePeakDist, the peak will not
81  * be added and the closest peak will be flagged as detected by the filter defined in keys.
82  *
83  * If the footprint does not overlap it will do nothing.
84  */
85  void add(std::shared_ptr<Footprint> footprint, afw::table::SchemaMapper const &peakSchemaMapper,
86  KeyTuple const &keys, float minNewPeakDist = -1., float maxSamePeakDist = -1.) {
87  if (addSpans(footprint)) {
88  _footprints.push_back(footprint);
89  _source->set(keys.footprint, true);
90  _addPeaks(footprint->getPeaks(), &peakSchemaMapper, &keys, minNewPeakDist, maxSamePeakDist, nullptr);
91  }
92  }
93 
94  /*
95  * Merge an already-merged clump of Footprints into this
96  *
97  * If minNewPeakDist >= 0, it will add all peaks from the footprint to the merged Footprint
98  * that are greater than minNewPeakDist away from the closest peak in the existing list.
99  * If minNewPeakDist < 0, no peaks will be added from the footprint.
100  *
101  * If maxSamePeakDist >= 0, it will find the closest peak in the existing list for every peak
102  * in the footprint. If the closest peak is less than maxSamePeakDist, the peak will not
103  * be added to the list and the flags from the closest peak will be set to the OR of the two.
104  *
105  * If the FootprintMerge does not overlap it will do nothing.
106  */
107  void add(FootprintMerge const &other, FilterMap const &keys, float minNewPeakDist = -1.,
108  float maxSamePeakDist = -1.) {
109  if (addSpans(other.getMergedFootprint())) {
110  _footprints.insert(_footprints.end(), other._footprints.begin(), other._footprints.end());
111  // Set source flags to the OR of the flags of the two inputs
112  for (auto const &key : keys) {
113  afw::table::Key<afw::table::Flag> const &flagKey = key.second.footprint;
114  _source->set(flagKey, _source->get(flagKey) || other._source->get(flagKey));
115  }
116  _addPeaks(other.getMergedFootprint()->getPeaks(), nullptr, nullptr, minNewPeakDist, maxSamePeakDist,
117  &keys);
118  }
119  }
120 
121  // Get the bounding box of the merge
122  lsst::geom::Box2I getBBox() const { return getMergedFootprint()->getBBox(); }
123 
124  std::shared_ptr<Footprint> getMergedFootprint() const { return _source->getFootprint(); }
125 
127 
128  // Implementation helper for add() methods; returns true if the Footprint actually overlapped
129  // and was merged, and false otherwise.
131  if (!getMergedFootprint()->getSpans()->overlaps(*(footprint->getSpans()))) return false;
132  getMergedFootprint()->setSpans(getMergedFootprint()->getSpans()->union_(*(footprint->getSpans())));
133  return true;
134  }
135 
136 private:
137  /*
138  * Add new peaks to the list of peaks of the merged footprint.
139  * This function handles two different cases:
140  * - The peaks come from a single footprint. In this case, the peakSchemaMapper
141  * and keys should be defined so that it can create a new peak, copy the appropriate
142  * data, and set the peak flag defined in keys.
143  * - The peaks come from another FootprintMerge. In this case, filterMap should
144  * be defined so that the information from the other peaks can be propagated.
145  */
146  void _addPeaks(PeakCatalog const &otherPeaks, afw::table::SchemaMapper const *peakSchemaMapper,
147  KeyTuple const *keys, float minNewPeakDist, float maxSamePeakDist,
148  FilterMap const *filterMap) {
149  if (minNewPeakDist < 0 && maxSamePeakDist < 0) return;
150 
151  assert(peakSchemaMapper || filterMap);
152 
153  PeakCatalog &currentPeaks = getMergedFootprint()->getPeaks();
154  std::shared_ptr<PeakRecord> nearestPeak;
155  // Create new list of peaks
156  PeakCatalog newPeaks(currentPeaks.getTable());
157  float minNewPeakDist2 = minNewPeakDist * minNewPeakDist;
158  float maxSamePeakDist2 = maxSamePeakDist * maxSamePeakDist;
159  for (PeakCatalog::const_iterator otherIter = otherPeaks.begin(); otherIter != otherPeaks.end();
160  ++otherIter) {
161  float minDist2 = std::numeric_limits<float>::infinity();
162 
163  for (PeakCatalog::const_iterator currentIter = currentPeaks.begin();
164  currentIter != currentPeaks.end(); ++currentIter) {
165  float dist2 = otherIter->getI().distanceSquared(currentIter->getI());
166 
167  if (dist2 < minDist2) {
168  minDist2 = dist2;
169  nearestPeak = currentIter;
170  }
171  }
172 
173  if (minDist2 < maxSamePeakDist2 && nearestPeak && maxSamePeakDist > 0) {
174  if (peakSchemaMapper) {
175  nearestPeak->set(keys->peak, true);
176  } else {
177  for (auto const &i : *filterMap) {
178  afw::table::Key<afw::table::Flag> const &flagKey = i.second.peak;
179  nearestPeak->set(flagKey, nearestPeak->get(flagKey) || otherIter->get(flagKey));
180  }
181  }
182  } else if (minDist2 > minNewPeakDist2 && !(minNewPeakDist < 0)) {
183  if (peakSchemaMapper) {
184  std::shared_ptr<PeakRecord> newPeak = newPeaks.addNew();
185  newPeak->assign(*otherIter, *peakSchemaMapper);
186  newPeak->set(keys->peak, true);
187  } else {
188  newPeaks.push_back(otherIter);
189  }
190  }
191  }
192 
193  getMergedFootprint()->getPeaks().insert(getMergedFootprint()->getPeaks().end(), newPeaks.begin(),
194  newPeaks.end(),
195  true // deep-copy
196  );
197  }
198 
201 };
202 
204  std::vector<std::string> const &filterList,
205  afw::table::Schema const &initialPeakSchema)
206  : _peakSchemaMapper(initialPeakSchema) {
207  _initialize(sourceSchema, filterList);
208 }
209 
211  std::vector<std::string> const &filterList)
212  : _peakSchemaMapper(PeakTable::makeMinimalSchema()) {
213  _initialize(sourceSchema, filterList);
214 }
215 
221 
222 void FootprintMergeList::_initialize(afw::table::Schema &sourceSchema,
223  std::vector<std::string> const &filterList) {
224  _peakSchemaMapper.addMinimalSchema(_peakSchemaMapper.getInputSchema(), true);
225  // Add Flags for the filters
226  for (auto const &iter : filterList) {
227  KeyTuple &keys = _filterMap[iter];
228  keys.footprint = sourceSchema.addField<afw::table::Flag>(
229  "merge_footprint_" + iter,
230  "Detection footprint overlapped with a detection from filter " + iter);
231  keys.peak = _peakSchemaMapper.editOutputSchema().addField<afw::table::Flag>(
232  "merge_peak_" + iter, "Peak detected in filter " + iter);
233  }
234  _peakTable = PeakTable::make(_peakSchemaMapper.getOutputSchema());
235 }
236 
238  afw::table::SourceCatalog const &inputCat, std::string const &filter,
239  float minNewPeakDist, bool doMerge, float maxSamePeakDist) {
240  FilterMap::const_iterator keyIter = _filterMap.find(filter);
241  if (keyIter == _filterMap.end()) {
243  (boost::format("Filter %s not in original list") % filter).str());
244  }
245 
246  // If list is empty or merging not requested, don't check for any matches, just add all the objects
247  bool checkForMatches = !_mergeList.empty() && doMerge;
248 
249  for (afw::table::SourceCatalog::const_iterator srcIter = inputCat.begin(); srcIter != inputCat.end();
250  ++srcIter) {
251  // Only consider unblended objects
252  if (srcIter->getParent() != 0) continue;
253 
254  std::shared_ptr<Footprint> foot = srcIter->getFootprint();
255 
256  // Empty pointer to account for the first match in the catalog. If there is more than one
257  // match, subsequent matches will be merged with this one
259 
260  if (checkForMatches) {
261  FootprintMergeVec::iterator iter = _mergeList.begin();
262  while (iter != _mergeList.end()) {
263  // Grow by one pixel to allow for touching
264  lsst::geom::Box2I box((**iter).getBBox());
265  box.grow(lsst::geom::Extent2I(1, 1));
266  if (box.overlaps(foot->getBBox()) && (**iter).overlaps(*foot)) {
267  if (!first) {
268  first = *iter;
269  // Spatially extend existing FootprintMerge in order to connect subsequent,
270  // now-overlapping FootprintMerges. If a subsequent FootprintMerge overlaps with
271  // the new footprint, it's now guaranteed to overlap with this first FootprintMerge.
272  // Hold off adding foot's lower-priority footprints and peaks until the
273  // higher-priority existing peaks are merged into this first FootprintMerge.
274  first->addSpans(foot);
275  } else {
276  // Add existing merged Footprint to first
277  first->add(**iter, _filterMap, minNewPeakDist, maxSamePeakDist);
278  iter = _mergeList.erase(iter);
279  continue;
280  }
281  }
282  ++iter;
283  } // while mergeList
284  } // if checkForMatches
285 
286  if (first) {
287  // Now merge footprint including peaks into the newly-connected, higher-priority FootprintMerge
288  first->add(foot, _peakSchemaMapper, keyIter->second, minNewPeakDist, maxSamePeakDist);
289  } else {
290  // Footprint did not overlap with any existing FootprintMerges. Add to MergeList
291  _mergeList.push_back(std::make_shared<FootprintMerge>(foot, sourceTable, _peakTable,
292  _peakSchemaMapper, keyIter->second));
293  }
294  }
295 }
296 
298  // Now set the merged footprint as the footprint of the SourceRecord
299  for (auto const &iter : _mergeList) {
300  outputCat.push_back((*iter).getSource());
301  }
302 }
303 } // namespace detection
304 } // namespace afw
305 } // namespace lsst
int end
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition: Exception.h:48
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)
FootprintMergeList::FilterMap FilterMap
std::shared_ptr< afw::table::SourceRecord > getSource() const
std::shared_ptr< Footprint > getMergedFootprint() const
FootprintMerge(FootprintMerge &&)=default
FootprintMergeList::KeyTuple KeyTuple
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.)
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:96
Iterator class for CatalogT.
Definition: Catalog.h:40
iterator begin()
Iterator access.
Definition: Catalog.h:401
std::shared_ptr< Table > getTable() const
Return the table associated with the catalog.
Definition: Catalog.h:115
Defines the fields and offsets for a table.
Definition: Schema.h:51
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:479
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
typename Base::const_iterator const_iterator
Definition: SortedCatalog.h:50
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)