Belle II Software  release-05-02-19
CDCSensitiveDetector.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Guofu Cao *
7  * 2012.03.05 SensitiveDetector -> CDCSensitiveDetector by M. Uchida *
8  * *
9  * This software is provided "as is" without any warranty. *
10  **************************************************************************/
11 
12 #ifndef CDCSensitiveDetector_H
13 #define CDCSensitiveDetector_H
14 
15 #include <simulation/kernel/SensitiveDetectorBase.h>
16 #include <cdc/dataobjects/WireID.h>
17 
18 #include <vector>
19 #include <map>
20 
21 
22 namespace Belle2 {
27  class CDCSimHit;
28  namespace CDC {
29 
30  class CDCGeometryPar;
31 
33 
37  class CDCSensitiveDetector: public Simulation::SensitiveDetectorBase {
38 
39  public:
40 
42  CDCSensitiveDetector(G4String name, G4double thresholdEnergyDeposit, G4double thresholdKineticEnergy);
43 
46 
48  void Initialize(G4HCofThisEvent*) override;
49 
51  bool step(G4Step* aStep, G4TouchableHistory* history) override;
52 
54  // void BeginOfEvent(G4HCofThisEvent*);
55 
57  void EndOfEvent(G4HCofThisEvent*) override;
58 
60  void saveSimHit(const G4int layerId,
61  const G4int wireId,
62  const G4int trackID,
63  const G4int pid,
64  const G4double distance,
65  const G4double tof,
66  const G4double edep,
67  const G4double stepLength,
68  const G4ThreeVector& mom,
69  const G4ThreeVector& posW,
70  const G4ThreeVector& posIn,
71  const G4ThreeVector& posOut,
72  const G4ThreeVector& posTrack,
73  const G4int lr,
74  const G4int NewLrRaw,
75  const G4int NewLr,
76  const G4double speed,
77  const G4double hitWeight);
78 
79  //void AddbgOne(bool doit);
80 
81  private:
82 
84  G4double ClosestApproach(G4ThreeVector bwp, G4ThreeVector fwp, G4ThreeVector posIn, G4ThreeVector posOut,
85  G4ThreeVector& hitPosition, G4ThreeVector& wirePosition);//,G4double& transferT);
86 
88 
107  void CellBound(const G4int layerId, const G4int ic1, const G4int ic2,
108  const G4double venter[6], const G4double vexit[6],
109  const G4double s1, const G4double s2,
110  G4double xint[6], G4double& sint, G4int& iflag);
111 
112 
114 
128  void GCUBS(const G4double x, const G4double y, const G4double d1, const G4double d2, G4double a[4]);
129 
137  void for_Rotat(const G4double bfld[3]);
138 
145  void Rotat(G4double& x, G4double& y,
146  G4double& z,
147  const int mode);
148 
152  void Rotat(G4double x[3], const int mode);
153 
155 
180  void HELWIR(const G4double xwb4,
181  const G4double ywb4,
182  const G4double zwb4,
183  const G4double xwf4,
184  const G4double ywf4,
185  const G4double zwf4,
186  const G4double xp,
187  const G4double yp,
188  const G4double zp,
189  const G4double px,
190  const G4double py,
191  const G4double pz,
192  const G4double B_kG[3],
193  const G4double charge,
194  const G4int ntryMax,
195  G4double& distance,
196  G4double q2[3],
197  G4double q1[3],
198  G4double q3[3],
199  G4int& ntry);
200 
202 
214  void Mvopr(const G4int ndim, const G4double b[3], const G4double m[3][3],
215  const G4double a[3], G4double c[3], const G4int mode);
216 
218  std::vector<int> WireId_in_hit_order(int id0, int id1, int nWires);
219 
222 
224  void reAssignLeftRightInfo();
225 
231  unsigned short areNeighbors(const WireID& wireId, const WireID& otherWireId) const;
232 
241  unsigned short areNeighbors(unsigned short iCLayer, unsigned short iSuperLayer, unsigned short iLayer, unsigned short iWire,
242  const WireID& otherWireId) const;
243 
244 
248  G4int m_magneticField;
249 
254  G4int m_nonUniformField;
255 
256  G4double m_brot[3][3];
258  // CDCGeometryPar& m_cdcgp; /**< Reference to CDCGeometryPar object. */
264  G4double m_thresholdEnergyDeposit;
265 
269  G4double m_thresholdKineticEnergy;
270 
271  G4bool m_wireSag;
275  G4double m_minTrackLength;
277  int m_hitNumber;
279  std::multimap<unsigned short, CDCSimHit*>
282  std::vector<CDCSimHit*> m_hitWithNegWeight;
284  // std::vector<std::multimap<unsigned short, CDCSimHit*>::iterator> m_posWeightMapItBegin; /**< Iterator showing begin of positive hit map*/
285 
286  // std::vector<std::multimap<unsigned short, CDCSimHit*>::iterator> m_posWeightMapItEnd; /**< Iterator showing end of positive hit map*/
287 
288  // int m_nPosHits, m_nNegHits;
289 
290  const signed short CCW = 1;
291  const signed short CW = -1;
292  const signed short CW_OUT_NEIGHBOR = 1;
293  const signed short CW_NEIGHBOR = 3;
294  const signed short CW_IN_NEIGHBOR = 5;
295  const signed short CCW_IN_NEIGHBOR = 7;
296  const signed short CCW_NEIGHBOR = 9;
297  const signed short CCW_OUT_NEIGHBOR = 11;
298  };
299  }
301 } // end of namespace Belle2
302 
303 #endif
Belle2::WireID
Class to identify a wire inside the CDC.
Definition: WireID.h:44
Belle2::CDC::CDCSensitiveDetector::Mvopr
void Mvopr(const G4int ndim, const G4double b[3], const G4double m[3][3], const G4double a[3], G4double c[3], const G4int mode)
Calculate the result of a matrix times vector.
Definition: CDCSensitiveDetector.cc:1339
Belle2::CDC::CDCSensitiveDetector::CW_IN_NEIGHBOR
const signed short CW_IN_NEIGHBOR
Constant for clockwise inwards.
Definition: CDCSensitiveDetector.h:303
Belle2::CDC::CDCSensitiveDetector::WireId_in_hit_order
std::vector< int > WireId_in_hit_order(int id0, int id1, int nWires)
Sort wire id.
Definition: CDCSensitiveDetector.cc:1388
Belle2::CDC::CDCSensitiveDetector::CW
const signed short CW
Constant for clockwise orientation.
Definition: CDCSensitiveDetector.h:300
Belle2::CDC::CDCSensitiveDetector::m_thresholdKineticEnergy
G4double m_thresholdKineticEnergy
Threshold kinetic energy to be stored.
Definition: CDCSensitiveDetector.h:278
Belle2::CDC::CDCSensitiveDetector::m_magneticField
G4int m_magneticField
Magnetic field is on or off.
Definition: CDCSensitiveDetector.h:257
Belle2::CDC::CDCSensitiveDetector::CellBound
void CellBound(const G4int layerId, const G4int ic1, const G4int ic2, const G4double venter[6], const G4double vexit[6], const G4double s1, const G4double s2, G4double xint[6], G4double &sint, G4int &iflag)
Calculate intersection of track with cell boundary.
Definition: CDCSensitiveDetector.cc:696
Belle2::CDC::CDCSensitiveDetector::CCW
const signed short CCW
Constant for counterclockwise orientation.
Definition: CDCSensitiveDetector.h:299
Belle2::CDC::CDCSensitiveDetector::GCUBS
void GCUBS(const G4double x, const G4double y, const G4double d1, const G4double d2, G4double a[4])
Definition: CDCSensitiveDetector.cc:943
Belle2::CDC::CDCSensitiveDetector::reAssignLeftRightInfo
void reAssignLeftRightInfo()
Re-assign left/right info.
Definition: CDCSensitiveDetector.cc:1532
Belle2::CDC::CDCSensitiveDetector::CDCSensitiveDetector
CDCSensitiveDetector(G4String name, G4double thresholdEnergyDeposit, G4double thresholdKineticEnergy)
Constructor.
Definition: CDCSensitiveDetector.cc:53
Belle2::CDC::CDCSensitiveDetector::Rotat
void Rotat(G4double &x, G4double &y, G4double &z, const int mode)
Translation method.
Definition: CDCSensitiveDetector.cc:1014
Belle2::CDC::CDCSensitiveDetector::step
bool step(G4Step *aStep, G4TouchableHistory *history) override
Process each step and calculate variables defined in CDCB4VHit.
Definition: CDCSensitiveDetector.cc:122
Belle2::CDC::CDCSensitiveDetector::CW_OUT_NEIGHBOR
const signed short CW_OUT_NEIGHBOR
Constant for clockwise outwards.
Definition: CDCSensitiveDetector.h:301
Belle2::CDC::CDCSensitiveDetector::m_hitWithNegWeight
std::vector< CDCSimHit * > m_hitWithNegWeight
Vector containing hits with negative weight.
Definition: CDCSensitiveDetector.h:291
Belle2::CDC::CDCGeometryPar
The Class for CDC Geometry Parameters.
Definition: CDCGeometryPar.h:75
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::CDC::CDCSensitiveDetector::m_hitNumber
int m_hitNumber
The current number of created hits in an event.
Definition: CDCSensitiveDetector.h:286
Belle2::CDC::CDCSensitiveDetector::m_cdcgp
CDCGeometryPar * m_cdcgp
Pointer to CDCGeometryPar object.
Definition: CDCSensitiveDetector.h:268
Belle2::CDC::CDCSensitiveDetector::for_Rotat
void for_Rotat(const G4double bfld[3])
Calculates a rotation matrix.
Definition: CDCSensitiveDetector.cc:978
Belle2::CDC::CDCSensitiveDetector::saveSimHit
void saveSimHit(const G4int layerId, const G4int wireId, const G4int trackID, const G4int pid, const G4double distance, const G4double tof, const G4double edep, const G4double stepLength, const G4ThreeVector &mom, const G4ThreeVector &posW, const G4ThreeVector &posIn, const G4ThreeVector &posOut, const G4ThreeVector &posTrack, const G4int lr, const G4int NewLrRaw, const G4int NewLr, const G4double speed, const G4double hitWeight)
Save CDCSimHit into datastore.
Definition: CDCSensitiveDetector.cc:589
Belle2::CDC::CDCSensitiveDetector::CCW_OUT_NEIGHBOR
const signed short CCW_OUT_NEIGHBOR
Constant for counterclockwise outwards.
Definition: CDCSensitiveDetector.h:306
Belle2::CDC::CDCSensitiveDetector::m_nonUniformField
G4int m_nonUniformField
Magnetic field is uniform or non-uniform.
Definition: CDCSensitiveDetector.h:263
Belle2::CDC::CDCSensitiveDetector::m_thresholdEnergyDeposit
G4double m_thresholdEnergyDeposit
Threshold energy deposit to be stored.
Definition: CDCSensitiveDetector.h:273
Belle2::CDC::CDCSensitiveDetector::m_minTrackLength
G4double m_minTrackLength
Min.
Definition: CDCSensitiveDetector.h:284
Belle2::CDC::CDCSensitiveDetector::CCW_NEIGHBOR
const signed short CCW_NEIGHBOR
Constant for counterclockwise.
Definition: CDCSensitiveDetector.h:305
Belle2::CDC::CDCSensitiveDetector::setModifiedLeftRightFlag
void setModifiedLeftRightFlag()
set left/right flag modified for tracking
Definition: CDCSensitiveDetector.cc:1446
Belle2::CDC::CDCSensitiveDetector::HELWIR
void HELWIR(const G4double xwb4, const G4double ywb4, const G4double zwb4, const G4double xwf4, const G4double ywf4, const G4double zwf4, const G4double xp, const G4double yp, const G4double zp, const G4double px, const G4double py, const G4double pz, const G4double B_kG[3], const G4double charge, const G4int ntryMax, G4double &distance, G4double q2[3], G4double q1[3], G4double q3[3], G4int &ntry)
Calculate closest points between helix and wire.
Definition: CDCSensitiveDetector.cc:1069
Belle2::CDC::CDCSensitiveDetector::CCW_IN_NEIGHBOR
const signed short CCW_IN_NEIGHBOR
Constant for counterclockwise inwards.
Definition: CDCSensitiveDetector.h:304
Belle2::CDC::CDCSensitiveDetector::Initialize
void Initialize(G4HCofThisEvent *) override
Register CDC hits collection into G4HCofThisEvent.
Definition: CDCSensitiveDetector.cc:96
Belle2::CDC::CDCSensitiveDetector::areNeighbors
unsigned short areNeighbors(const WireID &wireId, const WireID &otherWireId) const
Check if neighboring cell in the same super-layer; essentially a copy from cdcLocalTracking/mclookup.
Definition: CDCSensitiveDetector.cc:1618
Belle2::CDC::CDCSensitiveDetector::EndOfEvent
void EndOfEvent(G4HCofThisEvent *) override
Do what you want to do at the beginning of each event (why this is not called ?)
Definition: CDCSensitiveDetector.cc:583
Belle2::CDC::CDCSensitiveDetector::m_hitWithPosWeight
std::multimap< unsigned short, CDCSimHit * > m_hitWithPosWeight
Map containing hits with positive weight.
Definition: CDCSensitiveDetector.h:289
Belle2::CDC::CDCSensitiveDetector::~CDCSensitiveDetector
~CDCSensitiveDetector()
Destructor.
Definition: CDCSensitiveDetector.h:54
Belle2::CDC::CDCSensitiveDetector::ClosestApproach
G4double ClosestApproach(G4ThreeVector bwp, G4ThreeVector fwp, G4ThreeVector posIn, G4ThreeVector posOut, G4ThreeVector &hitPosition, G4ThreeVector &wirePosition)
Assume line track to calculate distance between track and wire (drift length).
Definition: CDCSensitiveDetector.cc:1419
Belle2::CDC::CDCSensitiveDetector::m_wireSag
G4bool m_wireSag
Switch to activate wire sag effect.
Definition: CDCSensitiveDetector.h:280
Belle2::CDC::CDCSensitiveDetector::CW_NEIGHBOR
const signed short CW_NEIGHBOR
Constant for clockwise.
Definition: CDCSensitiveDetector.h:302
Belle2::CDC::CDCSensitiveDetector::m_modifiedLeftRightFlag
G4bool m_modifiedLeftRightFlag
Switch for left/right flag modified for tracking.
Definition: CDCSensitiveDetector.h:282
Belle2::CDC::CDCSensitiveDetector::m_brot
G4double m_brot[3][3]
a rotation matrix.
Definition: CDCSensitiveDetector.h:265