Belle II Software development
ReattachCDCWireHitsToRecoTracksModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8#include <tracking/modules/reattachCDCWireHitsToRecoTracks/ReattachCDCWireHitsToRecoTracksModule.h>
9
10#include <tracking/trackFindingCDC/topology/CDCWire.h>
11#include <tracking/trackFindingCDC/geometry/Vector3D.h>
12#include <tracking/trackFindingCDC/geometry/Vector2D.h>
13#include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectory3D.h>
14#include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectory2D.h>
15#include <tracking/trackFindingCDC/eventdata/trajectories/CDCTrajectorySZ.h>
16#include <tracking/trackFitting/fitter/base/TrackFitter.h>
17#include <tracking/dataobjects/RecoHitInformation.h>
18#include <tracking/dbobjects/DAFConfiguration.h>
19
20using namespace Belle2;
21using namespace TrackFindingCDC;
22
23REG_MODULE(ReattachCDCWireHitsToRecoTracks);
24
26 Module()
27{
28 setDescription(R"DOC(
29Module to loop over low-ADC/TOT CDCWireHits and RecoTracks
30and reattach the hits to the tracks if they are closer
31than a given distance.)DOC");
32 setPropertyFlags(c_ParallelProcessingCertified); // parallel processing
33
34 addParam("CDCWireHitsStoreArrayName", m_CDCWireHitsStoreArrayName,
35 "Name of the input CDCWireHit StoreArray", m_CDCWireHitsStoreArrayName);
36 addParam("inputRecoTracksStoreArrayName", m_inputRecoTracksStoreArrayName,
37 "Name of the input RecoTrack StoreArray", m_inputRecoTracksStoreArrayName);
38 addParam("outputRecoTracksStoreArrayName", m_outputRecoTracksStoreArrayName,
39 "Name of the output RecoTrack StoreArray", m_outputRecoTracksStoreArrayName);
40 addParam("MaximumDistance", m_maximumDistance,
41 "Distance (cm) below which (exclusive) a CDC hit can be reattached to a track", m_maximumDistance);
42 addParam("MinimumADC", m_minimumADC,
43 "ADC above which (inclusive) a CDC hit can be reattached to a track", m_minimumADC);
44 addParam("MinimumTOT", m_minimumTOT,
45 "TOT above which (inclusive) a CDC hit can be reattached to a track", m_minimumTOT);
46 addParam("MaximumAbsD0", m_maximumAbsD0,
47 "Only tracks with an absolute value of d0 below (exclusive) this parameter (cm) are considered", m_maximumAbsD0);
48 addParam("MaximumAbsZ0", m_maximumAbsZ0,
49 "Only tracks with an absolute value of z0 below (exclusive) this parameter (cm) are considered", m_maximumAbsZ0);
50 addParam("trackFitType", m_trackFitType,
51 "Type of track fit algorithm to use the corresponding DAFParameter, the list is defined in DAFConfiguration class.",
53
54}
55
56
64
65
73
75{
76
78
79 for (RecoTrack& recoTrack : m_inputRecoTracks) {
80 // only fit tracks coming from the IP (d0 and z0 from Helix)
81 const Vector3D trackPosition(recoTrack.getPositionSeed());
82 const Vector3D trackMomentum(recoTrack.getMomentumSeed());
83 const CDCTrajectory3D trajectory(trackPosition, recoTrack.getTimeSeed(), trackMomentum, recoTrack.getChargeSeed());
84 const CDCTrajectory2D& trajectory2D(trajectory.getTrajectory2D());
85 const CDCTrajectorySZ& trajectorySZ(trajectory.getTrajectorySZ());
86 const double d0Estimate((trajectory2D.getClosest(Vector2D(0, 0))).norm());
87 const double z0Estimate(trajectorySZ.getZ0());
88 if (abs(d0Estimate) < m_maximumAbsD0 and abs(z0Estimate) < m_maximumAbsZ0) {
89 if (trackFitter.fit(recoTrack)) {
90 m_mapToHitsOnTrack[&recoTrack] = recoTrack.getSortedCDCHitList();
91 }
92 }
93 }
94
95 // Loop over the CDC hits and find the closest track (if any) whose distance to the hit is smaller than the threshold.
96 // Only the hits with the BadADCOrTOTFlag are considered (these are hits rejected by the TFCDC_WireHitPreparer module).
97 for (CDCWireHit& wireHit : *m_CDCWireHits) {
98 if ((wireHit.getAutomatonCell().hasBadADCOrTOTFlag()) and
99 (wireHit.getHit()->getADCCount() >= m_minimumADC) and
100 (wireHit.getHit()->getTOT() >= m_minimumTOT)) {
101
102 double currentMinimalDistance(m_maximumDistance);
103 RecoTrack* currentClosestTrack(nullptr);
104 ERightLeft currentRlInfo(ERightLeft::c_Unknown);
105
106 for (RecoTrack& recoTrack : m_inputRecoTracks) {
107 if (m_mapToHitsOnTrack.find(&recoTrack) == m_mapToHitsOnTrack.end()) { // Track not considered
108 continue;
109 }
110
111 bool neighborFound(false);
112 for (CDCHit* hitOnTrack : m_mapToHitsOnTrack[&recoTrack]) {
113
114 if (neighborFound) {
115 continue;
116 }
117 // To be added, the hit off track needs to be a neighbor of at least one hit on track
118 if (wireHit.getWire().isPrimaryNeighborWith(*CDCWire::getInstance(*hitOnTrack))) {
119 neighborFound = true;
120 } else {
121 continue;
122 }
123
124 const ReconstructionResults results(reconstruct(wireHit,
125 recoTrack,
126 recoTrack.getRecoHitInformation(hitOnTrack)));
127
128 if (not results.isValid) {
129 continue;
130 }
131
132 if (std::abs(results.distanceToTrack) < currentMinimalDistance) {
133
134 currentMinimalDistance = std::abs(results.distanceToTrack);
135 currentClosestTrack = &recoTrack;
136 currentRlInfo = results.rlInfo;
137
138 B2DEBUG(29, "Background hit close to the track found..." << std::endl
139 << "Layer of the hit on track: " << hitOnTrack->getICLayer() << std::endl
140 << "Layer of the background hit: " << wireHit.getHit()->getICLayer() << std::endl
141 << "ID of the background hit: " << wireHit.getHit()->getID() << std::endl
142 << "ADC of the background hit: " << wireHit.getHit()->getADCCount() << std::endl
143 << "TOT of the background hit: " << wireHit.getHit()->getTOT() << std::endl
144 << "Distance from track to hit: " << results.distanceToTrack << std::endl);
145 }
146 }
147 }
148
149 if (currentMinimalDistance < m_maximumDistance) { // This hit needs to be added to a RecoTrack
150 HitToAddInfo hitToAddInfo;
151 hitToAddInfo.hit = &wireHit;
152 hitToAddInfo.rlInfo = currentRlInfo;
153 m_mapToHitsToAdd[currentClosestTrack].emplace_back(hitToAddInfo);
154 }
155 }
156 }
157}
158
159
161{
162
163 for (RecoTrack& recoTrack : m_inputRecoTracks) {
164
165 RecoTrack* newRecoTrack = recoTrack.copyToStoreArray(m_outputRecoTracks);
166
167 if (m_mapToHitsToAdd.find(&recoTrack) == m_mapToHitsToAdd.end()) { // No hit to add
168
169 newRecoTrack->addHitsFromRecoTrack(&recoTrack);
170
171 } else { // At least one hit to add
172
173 std::unordered_map<CDCWireHit*, double> previousArcLength;
174 std::unordered_map<CDCWireHit*, double> currentArcLength;
175 // Initialise the arc-length maps to zero and unset the taken and background flags.
176 for (HitToAddInfo& hitToAddInfo : m_mapToHitsToAdd[&recoTrack]) {
177 previousArcLength[hitToAddInfo.hit] = 0.0;
178 currentArcLength[hitToAddInfo.hit] = 0.0;
179 (hitToAddInfo.hit)->getAutomatonCell().setTakenFlag(false);
180 (hitToAddInfo.hit)->getAutomatonCell().setBackgroundFlag(false);
181 }
182
183 unsigned int sortingParameter(0);
184 for (CDCHit* hitOnTrack : m_mapToHitsOnTrack[&recoTrack]) {
185 for (HitToAddInfo& hitToAddInfo : m_mapToHitsToAdd[&recoTrack]) {
186 CDCWireHit& hitToAdd = *(hitToAddInfo.hit);
187 if (not hitToAdd.getAutomatonCell().hasTakenFlag()) {
188
189 const ReconstructionResults results(reconstruct(hitToAdd,
190 recoTrack,
191 recoTrack.getRecoHitInformation(hitOnTrack)));
192
193 previousArcLength[&hitToAdd] = currentArcLength[&hitToAdd];
194 currentArcLength[&hitToAdd] = results.arcLength;
195
196 B2DEBUG(29, "Hit to be added..." << std::endl
197 << "Layer of the hit on track: " << hitOnTrack->getICLayer() << std::endl
198 << "Layer of the background hit: " << hitToAdd.getHit()->getICLayer() << std::endl
199 << "ID of the background hit: " << hitToAdd.getHit()->getID() << std::endl
200 << "ADC of the background hit: " << hitToAdd.getHit()->getADCCount() << std::endl
201 << "TOT of the background hit: " << hitToAdd.getHit()->getTOT() << std::endl
202 << "Distance from track to hit: " << results.distanceToTrack << std::endl
203 << "Previous arc length of the hit: " << previousArcLength[&hitToAdd] << std::endl
204 << "Current arc length of the hit: " << currentArcLength[&hitToAdd] << std::endl);
205
206 if ((previousArcLength[&hitToAdd] > 0) and (currentArcLength[&hitToAdd] < 0)) { // Hit needs to be added here.
207
209 newRecoTrack->addCDCHit(hitToAdd.getHit(), sortingParameter, rl, RecoHitInformation::c_ReattachCDCWireHitsToRecoTracks);
210 hitToAdd.getAutomatonCell().setTakenFlag(true);
211 ++sortingParameter;
212
213 }
214 }
215 }
216 const RecoHitInformation::RightLeftInformation rl = recoTrack.getRecoHitInformation(hitOnTrack)->getRightLeftInformation();
217 const RecoHitInformation::OriginTrackFinder foundBy = recoTrack.getRecoHitInformation(hitOnTrack)->getFoundByTrackFinder();
218 newRecoTrack->addCDCHit(hitOnTrack, sortingParameter, rl, foundBy);
219 //TODO: In the (rare) case where more than one hit are added between the same 2 hits on track, one should order them w.r.t. the arcLength.
220 ++sortingParameter;
221 }
222 }
223 }
224}
225
226
228 const CDCWireHit& wireHit,
229 const RecoTrack& recoTrack, const RecoHitInformation* const recoHitInformation) const
230{
231 ReconstructionResults results;
232
233 try {
234
235 const genfit::MeasuredStateOnPlane& mSoP(recoTrack.getMeasuredStateOnPlaneFromRecoHit(recoHitInformation));
236 const Vector3D trackPosition(mSoP.getPos());
237 const Vector3D trackMomentum(mSoP.getMom());
238 const CDCTrajectory3D trajectory(trackPosition, mSoP.getTime(), trackMomentum, mSoP.getCharge());
239
240 const CDCTrajectory2D& trajectory2D(trajectory.getTrajectory2D());
241 const CDCTrajectorySZ& trajectorySZ(trajectory.getTrajectorySZ());
242
243 Vector2D recoPos2D;
244 if (wireHit.isAxial()) {
245 recoPos2D = wireHit.reconstruct2D(trajectory2D);
246 } else {
247 const CDCWire& wire(wireHit.getWire());
248 const Vector2D& posOnXYPlane(wireHit.reconstruct2D(trajectory2D));
249
250 const double arcLength(trajectory2D.calcArcLength2D(posOnXYPlane));
251 const double z(trajectorySZ.mapSToZ(arcLength));
252
253 const Vector2D& wirePos2DAtZ(wire.getWirePos2DAtZ(z));
254
255 const Vector2D& recoPosOnTrajectory(trajectory2D.getClosest(wirePos2DAtZ));
256 const double driftLength(wireHit.getRefDriftLength());
257 Vector2D disp2D(recoPosOnTrajectory - wirePos2DAtZ);
258 disp2D.normalizeTo(driftLength);
259 recoPos2D = wirePos2DAtZ + disp2D;
260 }
261
262
263 results.arcLength = trajectory2D.calcArcLength2D(recoPos2D);
264 results.z = trajectorySZ.mapSToZ(results.arcLength);
265 results.distanceToTrack = trajectory2D.getDist2D(recoPos2D);
266
267 const Vector3D hitPosition(wireHit.getWire().getWirePos3DAtZ(trackPosition.z()));
268
269 Vector3D trackPosToWire{hitPosition - trackPosition};
270 results.rlInfo = trackPosToWire.xy().isRightOrLeftOf(trackMomentum.xy());
271
272 results.isValid = true;
273
274 } catch (...) {
275 B2WARNING("Distance measurement does not work.");
276 results.isValid = false;
277 }
278
279 return results;
280}
281
282
284 ERightLeft rlInfo) const
285{
286 using RightLeftInformation = RecoHitInformation::RightLeftInformation;
287 if (rlInfo == ERightLeft::c_Left) {
288 return RightLeftInformation::c_left;
289 } else if (rlInfo == ERightLeft::c_Right) {
290 return RightLeftInformation::c_right;
291 } else if (rlInfo == ERightLeft::c_Invalid) {
292 return RightLeftInformation::c_invalidRightLeftInformation;
293 } else {
294 return RightLeftInformation::c_undefinedRightLeftInformation;
295 }
296}
Class containing the result of the unpacker in raw data and the result of the digitizer in simulation...
Definition CDCHit.h:40
unsigned short getICLayer() const
Getter for iCLayer (0-55).
Definition CDCHit.h:178
unsigned short getID() const
Getter for encoded wire number.
Definition CDCHit.h:193
unsigned short getADCCount() const
Getter for integrated charge.
Definition CDCHit.h:230
unsigned short getTOT() const
Getter for TOT.
Definition CDCHit.h:248
ETrackFitType
Enum for identifying the type of track fit algorythm ( or cosmic)
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
Module()
Constructor.
Definition Module.cc:30
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
This class stores additional information to every CDC/SVD/PXD hit stored in a RecoTrack.
OriginTrackFinder
The TrackFinder which has added the hit to the track.
RightLeftInformation
The RightLeft information of the hit which is only valid for CDC hits.
This is the Reconstruction Event-Data Model Track.
Definition RecoTrack.h:79
size_t addHitsFromRecoTrack(const RecoTrack *recoTrack, unsigned int sortingParameterOffset=0, bool reversed=false, std::optional< double > optionalMinimalWeight=std::nullopt)
Add all hits from another RecoTrack to this RecoTrack.
Definition RecoTrack.cc:240
bool addCDCHit(const UsedCDCHit *cdcHit, const unsigned int sortingParameter, RightLeftInformation rightLeftInformation=RightLeftInformation::c_undefinedRightLeftInformation, OriginTrackFinder foundByTrackFinder=OriginTrackFinder::c_undefinedTrackFinder)
Adds a cdc hit with the given information to the reco track.
Definition RecoTrack.h:243
RecoTrack * copyToStoreArray(StoreArray< RecoTrack > &storeArray) const
Append a new RecoTrack to the given store array and copy its general properties, but not the hits the...
Definition RecoTrack.cc:529
static void registerRequiredRelations(StoreArray< RecoTrack > &recoTracks, std::string const &pxdHitsStoreArrayName="", std::string const &svdHitsStoreArrayName="", std::string const &cdcHitsStoreArrayName="", std::string const &bklmHitsStoreArrayName="", std::string const &eklmHitsStoreArrayName="", std::string const &recoHitInformationStoreArrayName="")
Convenience method which registers all relations required to fully use a RecoTrack.
Definition RecoTrack.cc:53
const genfit::MeasuredStateOnPlane & getMeasuredStateOnPlaneFromRecoHit(const RecoHitInformation *recoHitInfo, const genfit::AbsTrackRep *representation=nullptr) const
Return genfit's MeasuredStateOnPlane on plane for associated with one RecoHitInformation.
Definition RecoTrack.cc:579
void setTakenFlag(bool setTo=true)
Sets the taken flag to the given value. Default value true.
bool hasTakenFlag() const
Gets the current state of the taken marker flag.
Particle trajectory as it is seen in xy projection represented as a circle.
double calcArcLength2D(const Vector2D &point) const
Calculate the travel distance from the start position of the trajectory.
Vector2D getClosest(const Vector2D &point) const
Calculates the closest approach on the trajectory to the given point.
double getDist2D(const Vector2D &point) const
Calculates the distance from the point to the trajectory as seen from the xy projection.
Particle full three dimensional trajectory.
CDCTrajectory2D getTrajectory2D() const
Getter for the two dimensional trajectory.
CDCTrajectorySZ getTrajectorySZ() const
Getter for the sz trajectory.
Linear trajectory in sz space.
double mapSToZ(const double s=0) const
Translates the travel distance to the z coordinate.
double getZ0() const
Getter for the z coordinate at zero travel distance.
Class representing a hit wire in the central drift chamber.
Definition CDCWireHit.h:55
const CDCHit * getHit() const
Getter for the CDCHit pointer into the StoreArray.
Definition CDCWireHit.h:159
double getRefDriftLength() const
Getter for the drift length at the reference position of the wire.
Definition CDCWireHit.h:224
bool isAxial() const
Indicator if the underlying wire is axial.
Definition CDCWireHit.h:197
const CDCWire & getWire() const
Getter for the CDCWire the hit is located on.
Definition CDCWireHit.h:168
AutomatonCell & getAutomatonCell() const
Mutable getter for the automaton cell.
Definition CDCWireHit.h:286
Vector2D reconstruct2D(const CDCTrajectory2D &trajectory2D) const
Reconstructs a position of primary ionisation on the drift circle.
Class representing a sense wire in the central drift chamber.
Definition CDCWire.h:58
Vector3D getWirePos3DAtZ(const double z) const
Gives position of the wire at the given z coordinate.
Definition CDCWire.h:196
Vector2D getWirePos2DAtZ(const double z) const
Gives the xy projected position of the wire at the given z coordinate.
Definition CDCWire.h:192
static const CDCWire * getInstance(const WireID &wireID)
Getter from the wireID convenience object. Does not construct a new object.
Definition CDCWire.cc:24
std::string m_outputRecoTracksStoreArrayName
Name of the output RecoTrack StoreArray.
std::string m_CDCWireHitsStoreArrayName
Name of the input CDCWireHit StoreWrappedObjPtr.
ReattachCDCWireHitsToRecoTracksModule()
Constructor of the module. Setting up parameters and description.
int m_minimumADC
ADC above which (inclusive) a CDC hit can be reattached to a track.
std::string m_inputRecoTracksStoreArrayName
Name of the input RecoTrack StoreArray.
std::unordered_map< RecoTrack *, std::vector< CDCHit * > > m_mapToHitsOnTrack
Map from a RecoTrack ptr to the vector of the hits that belong to this track.
RecoHitInformation::RightLeftInformation rightLeftInformationTranslator(ERightLeft rlInfo) const
Translate a TrackFindingCDC::ERightLeft into a RecoHitInformation::RightLeftInformation.
ReconstructionResults reconstruct(const CDCWireHit &wireHit, const RecoTrack &recoTrack, const RecoHitInformation *recoHitInformation) const
Compute distance from a CDCWireHit to a RecoTrack using the mSoP found with a RecoHitInformation.
double m_maximumDistance
Distance (cm) below which (exclusive) a CDC hit can be reattached to a track.
int m_minimumTOT
TOT above which (inclusive) a CDC hit can be reattached to a track.
StoreWrappedObjPtr< std::vector< CDCWireHit > > m_CDCWireHits
Input CDCWireHits.
double m_maximumAbsD0
Only tracks with an absolute value of d0 below (exclusive) this parameter (cm) are considered.
short m_trackFitType
Track Fit type to select the proper DAFParameter from DAFConfiguration; by default c_CDConly.
std::unordered_map< RecoTrack *, std::vector< HitToAddInfo > > m_mapToHitsToAdd
Map from a RecoTrack ptr to the vector of the hits that need to be added to this track.
double m_maximumAbsZ0
Only tracks with an absolute value of z0 below (exclusive) this parameter (cm) are considered.
A two dimensional vector which is equipped with functions for correct handling of orientation relate...
Definition Vector2D.h:32
double normalizeTo(const double toLength)
Normalizes the vector to the given length.
Definition Vector2D.h:313
A three dimensional vector.
Definition Vector3D.h:33
const Vector2D & xy() const
Getter for the xy projected vector ( reference ! )
Definition Vector3D.h:508
double z() const
Getter for the z coordinate.
Definition Vector3D.h:496
Algorithm class to handle the fitting of RecoTrack objects.
bool fit(RecoTrack &recoTrack, genfit::AbsTrackRep *trackRepresentation, bool resortHits=false) const
Fit a reco track with a given non-default track representation.
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Abstract base class for different kinds of events.
Internal structure to store the information about a hit to be added.