10 #include <tracking/trackFindingCDC/display/CDCSVGPlotter.h>
11 #include <tracking/trackFindingCDC/display/SVGPrimitivePlotter.h>
12 #include <tracking/trackFindingCDC/display/Styling.h>
14 #include <tracking/trackFindingCDC/rootification/StoreWrappedObjPtr.h>
16 #include <tracking/trackFindingCDC/filters/axialSegmentPair/MCAxialSegmentPairFilter.h>
17 #include <tracking/trackFindingCDC/filters/segmentPair/MCSegmentPairFilter.h>
18 #include <tracking/trackFindingCDC/filters/segmentTriple/MCSegmentTripleFilter.h>
19 #include <tracking/trackFindingCDC/mclookup/CDCMCSegment2DLookUp.h>
21 #include <tracking/trackFindingCDC/eventdata/tracks/CDCTrack.h>
22 #include <tracking/trackFindingCDC/eventdata/tracks/CDCSegmentTriple.h>
23 #include <tracking/trackFindingCDC/eventdata/tracks/CDCAxialSegmentPair.h>
24 #include <tracking/trackFindingCDC/eventdata/tracks/CDCSegmentPair.h>
25 #include <tracking/trackFindingCDC/eventdata/segments/CDCSegment2D.h>
26 #include <tracking/trackFindingCDC/eventdata/segments/CDCWireHitCluster.h>
27 #include <tracking/trackFindingCDC/eventdata/hits/CDCWireHit.h>
29 #include <tracking/trackFindingCDC/geometry/Vector3D.h>
30 #include <tracking/trackFindingCDC/geometry/Vector2D.h>
32 #include <tracking/trackFindingCDC/utilities/ReversedRange.h>
34 #include <tracking/dataobjects/RecoTrack.h>
36 #include <framework/datastore/StoreArray.h>
38 #include <cdc/dataobjects/CDCSimHit.h>
39 #include <cdc/dataobjects/CDCHit.h>
41 #include <mdst/dataobjects/MCParticle.h>
46 using namespace TrackFindingCDC;
49 template<
bool a_drawTrajectory,
class AObject>
50 struct DrawTrajectoryHelper;
52 template<
class AObject>
53 struct DrawTrajectoryHelper<true, AObject> {
55 static void draw(EventDataPlotter& plotter, AObject&
object,
const AttributeMap& attributeMap)
57 plotter.drawTrajectory(
object, attributeMap);
61 template<
class AObject>
62 struct DrawTrajectoryHelper<false, AObject> {
64 static void draw(EventDataPlotter& plotter, AObject&
object,
const AttributeMap& attributeMap)
66 plotter.draw(
object, attributeMap);
73 class FlightTimeOrder {
78 return (x->getFlightTime() < y->getFlightTime());
83 static void printDataStoreContent()
85 B2INFO(
"Current content of the DataStore:");
86 B2INFO(
"StoreArrays:");
91 B2INFO(
"StoreObjPtr:");
97 const AttributeMap c_defaultSVGAttributes({
99 {
"stroke-width",
"0.55"},
101 {
"transform",
"translate(0, 1120) scale(1,-1)"}
105 m_eventdataPlotter(std::make_unique<
SVGPrimitivePlotter>(c_defaultSVGAttributes), animate, forwardFade)
113 int default_width = 1120;
114 int default_height = 1120;
133 AttributeMap attributeMap;
134 attributeMap.emplace(
"stroke", stroke);
140 AttributeMap attributeMap;
141 attributeMap.emplace(
"stroke", stroke);
147 AttributeMap attributeMap;
148 attributeMap.emplace(
"stroke", stroke);
158 const std::string& stroke,
159 const std::string& strokeWidth)
162 if (stroke !=
"") styling.
setStroke(stroke);
164 drawStoreArray<const CDCHit>(storeArrayName, styling);
168 const std::string& stroke,
169 const std::string& strokeWidth)
171 B2INFO(
"Drawing simulated hits");
173 if (not storeArray) {
174 B2WARNING(
"StoreArray " << storeArrayName <<
" not present");
175 printDataStoreContent();
178 std::vector<CDCSimHit> simHitsRelatedToHits;
179 for (
const CDCHit& hit : storeArray) {
181 simHitsRelatedToHits.push_back(*hit.getRelated<
CDCSimHit>(
"CDCSimHits"));
183 B2INFO(
"#CDCSimHits: " << storeArray.getEntries());
195 const std::string& stroke,
196 const std::string& strokeWidth)
199 if (stroke !=
"") styling.
setStroke(stroke);
201 drawStoreVector<const CDCWireHitCluster>(storeObjName, styling);
205 const std::string& stroke,
206 const std::string& strokeWidth)
209 if (stroke !=
"") styling.
setStroke(stroke);
211 drawStoreVector<const CDCSegment2D>(storeObjName, styling);
215 const std::string& stroke,
216 const std::string& strokeWidth)
219 if (stroke !=
"") styling.
setStroke(stroke);
221 const bool drawTrajectories =
true;
222 drawStoreVector<const CDCSegment2D, drawTrajectories>(storeObjName, styling);
226 const std::string& stroke,
227 const std::string& strokeWidth)
230 if (stroke !=
"") styling.
setStroke(stroke);
232 drawStoreVector<const CDCAxialSegmentPair>(storeObjName, styling);
236 const std::string& stroke,
237 const std::string& strokeWidth)
240 if (stroke !=
"") styling.
setStroke(stroke);
242 drawStoreVector<const CDCSegmentPair>(storeObjName, styling);
246 const std::string& stroke,
247 const std::string& strokeWidth)
250 if (stroke !=
"") styling.
setStroke(stroke);
252 drawStoreVector<const CDCSegmentTriple>(storeObjName, styling);
256 const std::string& stroke,
257 const std::string& strokeWidth)
260 if (stroke !=
"") styling.
setStroke(stroke);
262 const bool drawTrajectories =
true;
263 drawStoreVector<const CDCSegmentTriple, drawTrajectories>(storeObjName, styling);
267 const std::string& stroke,
268 const std::string& strokeWidth)
271 if (stroke !=
"") styling.
setStroke(stroke);
273 drawStoreVector<const CDCTrack>(storeObjName, styling);
277 const std::string& stroke,
278 const std::string& strokeWidth)
281 if (stroke !=
"") styling.
setStroke(stroke);
283 const bool drawTrajectories =
true;
284 drawStoreVector<const CDCTrack, drawTrajectories>(storeObjName, styling);
288 const std::string& stroke,
289 const std::string& strokeWidth)
292 if (stroke !=
"") styling.
setStroke(stroke);
294 drawStoreArray<const RecoTrack>(storeArrayName, styling);
298 const std::string& stroke,
299 const std::string& strokeWidth)
302 if (stroke !=
"") styling.
setStroke(stroke);
304 const bool drawTrajectories =
true;
305 drawStoreArray<const RecoTrack, drawTrajectories>(storeArrayName, styling);
309 const std::string& stroke,
310 const std::string& strokeWidth)
313 if (stroke !=
"") styling.
setStroke(stroke);
315 const bool drawTrajectories =
true;
316 drawStoreArray<const MCParticle, drawTrajectories>(storeArrayName, styling);
320 const std::string& stroke,
321 const std::string& strokeWidth)
323 B2INFO(
"Drawing simulated hits connected by tof");
325 if (not hitStoreArray) {
326 B2WARNING(
"StoreArray " << hitStoreArrayName <<
" not present");
327 printDataStoreContent();
330 std::vector<CDCSimHit*> simHits;
331 for (
const CDCHit& hit : hitStoreArray) {
333 simHits.push_back(hit.getRelated<
CDCSimHit>());
337 std::map<int, std::set<CDCSimHit*, FlightTimeOrder>> simHitsByMcParticleId;
340 if (mcParticle !=
nullptr) {
342 simHitsByMcParticleId[mcTrackId].insert(simHit);
346 AttributeMap defaultAttributeMap = {{
"stroke", stroke}, {
"stroke-width", strokeWidth}};
348 for (
const auto& mcParticleIdAndSimHits : simHitsByMcParticleId) {
349 const std::set<CDCSimHit*, FlightTimeOrder>& simHitsForMcParticle =
350 mcParticleIdAndSimHits.second;
352 auto drawConnectSimHits = [
this, &defaultAttributeMap](
CDCSimHit * fromSimHit,
CDCSimHit * toSimHit) {
356 if (fromHit ==
nullptr)
return false;
357 if (toHit ==
nullptr)
return false;
366 Vector3D toDisplacement(toSimHit->getPosTrack() - toSimHit->getPosWire());
371 bool falseOrder =
false;
372 if (fromSimHit->
getArrayIndex() > toSimHit->getArrayIndex()) {
374 bool toReassigned = toHit->getRelatedWithWeight<
MCParticle>().second < 0;
375 if (not fromReassigned and not toReassigned) {
380 AttributeMap attributeMap = defaultAttributeMap;
382 attributeMap[
"stroke"] =
"red";
383 attributeMap[
"stroke-width"] =
"1.0";
385 draw(fromRecoHit2D, attributeMap);
386 draw(toRecoHit2D, attributeMap);
389 const float fromX = fromPos.
x();
390 const float fromY = fromPos.
y();
393 const float toX = toPos.
x();
394 const float toY = toPos.
y();
401 std::adjacent_find(simHitsForMcParticle.begin(),
402 simHitsForMcParticle.end(),
409 this->drawWrongRLHits<CDCSegment2D>(segmentsStoreObjName);
414 this->drawWrongRLHits<CDCTrack>(tracksStoreObjName);
417 template<
class ACDCHitCollection>
420 B2INFO(
"Draw wrong right left passage information from " << hitCollectionsStoreObjName);
422 if (not storedHitCollections) {
423 B2WARNING(hitCollectionsStoreObjName <<
"does not exist in current DataStore");
424 printDataStoreContent();
428 std::vector<ACDCHitCollection>& hitCollections = *storedHitCollections;
429 B2INFO(
"#HitCollections: " << hitCollections.size());
431 const CDCMCHitCollectionLookUp<ACDCHitCollection> mcHitCollectionLookUp;
434 for (
const ACDCHitCollection& hitCollection : hitCollections) {
435 EForwardBackward fbInfo = mcHitCollectionLookUp.isForwardOrBackwardToMCTrack(&hitCollection);
437 double rlPurity = mcHitCollectionLookUp.getRLPurity(&hitCollection);
438 int correctRLVote = mcHitCollectionLookUp.getCorrectRLVote(&hitCollection);
441 if (rlPurity < 0.5 and hitCollection.getAutomatonCell().hasAliasFlag())
continue;
444 if (correctRLVote < 0 and hitCollection.getAutomatonCell().hasReverseFlag())
continue;
447 for (
const auto& recoHit : hitCollection) {
448 ERightLeft rlInfo = recoHit.getRLInfo();
449 const CDCHit* hit = recoHit.getWireHit().getHit();
450 ERightLeft mcRLInfo = mcHitLookUp.getRLInfo(hit);
452 if (fbInfo == EForwardBackward::c_Backward) {
453 mcRLInfo = reversed(mcRLInfo);
456 std::string color =
"orange";
457 if (mcRLInfo != ERightLeft::c_Right and mcRLInfo != ERightLeft::c_Left) {
459 }
else if (mcRLInfo == rlInfo) {
461 }
else if (mcRLInfo == -rlInfo) {
465 AttributeMap attributeMap{{
"stroke", color}};
473 const std::string& stroke,
474 const std::string& strokeWidth)
476 B2INFO(
"Draw axial to axial segment pairs");
478 if (not storedSegments) {
479 B2WARNING(segmentsStoreObjName <<
"does not exist in current DataStore");
480 printDataStoreContent();
484 std::vector<CDCSegment2D>& segments = *storedSegments;
485 B2INFO(
"#Segments: " << segments.size());
487 std::vector<const CDCAxialSegment2D*> axialSegments;
489 if (segment.isAxial()) axialSegments.push_back(&segment);
493 std::vector<CDCAxialSegmentPair> mcAxialSegmentPairs;
496 if (fromSegment == toSegment)
continue;
497 mcAxialSegmentPairs.emplace_back(fromSegment, toSegment);
498 Weight mcWeight = mcAxialSegmentPairFilter(mcAxialSegmentPairs.back());
500 if (std::isnan(mcWeight)) {
501 mcAxialSegmentPairs.pop_back();
505 B2INFO(
"# Axial segment pairs: " << mcAxialSegmentPairs.size());
507 if (stroke !=
"") styling.
setStroke(stroke);
508 if (strokeWidth !=
"") styling.
setStroke(strokeWidth);
513 const std::string& stroke,
514 const std::string& strokeWidth)
516 B2INFO(
"Draw axial to stero segment pairs");
518 if (not storedSegments) {
519 B2WARNING(segmentsStoreObjName <<
"does not exist in current DataStore");
520 printDataStoreContent();
524 std::vector<CDCSegment2D>& segments = *storedSegments;
525 B2INFO(
"#Segments: " << segments.size());
527 std::vector<const CDCAxialSegment2D*> axialSegments;
528 std::vector<const CDCStereoSegment2D*> stereoSegments;
530 if (segment.isAxial()) {
531 axialSegments.push_back(&segment);
533 stereoSegments.push_back(&segment);
538 std::vector<CDCSegmentPair> mcSegmentPairs;
543 mcSegmentPairs.emplace_back(axialSegment, stereoSegment);
544 Weight mcWeight = mcSegmentPairFilter(mcSegmentPairs.back());
546 if (std::isnan(mcWeight)) {
547 mcSegmentPairs.pop_back();
552 mcSegmentPairs.emplace_back(stereoSegment, axialSegment);
553 Weight mcWeight = mcSegmentPairFilter(mcSegmentPairs.back());
555 if (std::isnan(mcWeight)) {
556 mcSegmentPairs.pop_back();
561 B2INFO(
"# Segment pairs: " << mcSegmentPairs.size());
563 if (stroke !=
"") styling.
setStroke(stroke);
564 if (strokeWidth !=
"") styling.
setStroke(strokeWidth);
569 const std::string& stroke,
570 const std::string& strokeWidth)
572 B2INFO(
"Draw segment triples");
574 if (not storedSegments) {
575 B2WARNING(segmentsStoreObjName <<
"does not exist in current DataStore");
576 printDataStoreContent();
580 std::vector<CDCSegment2D>& segments = *storedSegments;
581 B2INFO(
"#Segment " << segments.size());
583 std::vector<const CDCAxialSegment2D*> axialSegments;
584 std::vector<const CDCStereoSegment2D*> stereoSegments;
586 if (segment.isAxial()) {
587 axialSegments.push_back(&segment);
589 stereoSegments.push_back(&segment);
594 std::vector<CDCSegmentTriple> mcSegmentTriples;
598 if (startSegment == endSegment)
continue;
599 mcSegmentTriples.emplace_back(startSegment, middleSegment, endSegment);
600 Weight mcWeight = mcSegmentTripleFilter(mcSegmentTriples.back());
602 if (std::isnan(mcWeight)) {
603 mcSegmentTriples.pop_back();
608 B2INFO(
"# Segment triples: " << mcSegmentTriples.size());
610 if (stroke !=
"") styling.
setStroke(stroke);
611 if (strokeWidth !=
"") styling.
setStroke(strokeWidth);
620 float width = boundingBox.
getWidth();
622 int totalPoints = 1120 * 1120;
623 float svgHeight = roundf(sqrt(totalPoints * height / width));
624 float svgWidth = roundf(sqrt(totalPoints * width / height));
632 template <
class AItem,
bool a_drawTrajectories>
636 if (a_drawTrajectories) {
637 B2INFO(
"Drawing trajectories from StoreArray: " << storeArrayName);
639 B2INFO(
"Drawing StoreArray: " << storeArrayName);
642 using StoreItem =
typename std::remove_cv<AItem>::type;
644 if (not storeArray) {
645 B2WARNING(storeArrayName <<
" not present in the DataStore");
646 printDataStoreContent();
650 B2INFO(
"with " << storeArray.getEntries() <<
" entries");
651 drawIterable<a_drawTrajectories>(storeArray, styling);
652 B2INFO(
"Attributes are");
653 B2INFO(styling.
info());
656 template <
class AItem,
bool a_drawTrajectories>
658 Styling<AItem>& styling)
660 if (a_drawTrajectories) {
661 B2INFO(
"Drawing trajectories for vector from DataStore: " << storeObjName);
663 B2INFO(
"Drawing vector from DataStore: " << storeObjName);
666 using StoreItem =
typename std::remove_cv<AItem>::type;
668 if (not storeVector) {
669 B2WARNING(storeObjName <<
" not present in the DataStore");
670 B2INFO(
"Current content of the DataStore:");
671 printDataStoreContent();
675 const std::vector<StoreItem>& vector = *storeVector;
676 B2INFO(
"with " << vector.size() <<
" entries");
677 drawIterable<a_drawTrajectories>(reversedRange(vector), styling);
678 B2INFO(
"Attributes are");
679 B2INFO(styling.info());
682 template <
bool a_drawTrajectory,
class AIterable,
class AStyling>
685 unsigned int index = -1;
686 for (
const auto& item : items) {
689 AttributeMap attributeMap = styling.map(index, item);
690 draw<a_drawTrajectory>(item, attributeMap);
694 template <
bool a_drawTrajectory,
class AObject>