 |
Belle II Software
release-05-02-19
|
12 #ifndef CDCSensitiveDetector_H
13 #define CDCSensitiveDetector_H
15 #include <simulation/kernel/SensitiveDetectorBase.h>
16 #include <cdc/dataobjects/WireID.h>
37 class CDCSensitiveDetector:
public Simulation::SensitiveDetectorBase {
42 CDCSensitiveDetector(G4String name, G4double thresholdEnergyDeposit, G4double thresholdKineticEnergy);
51 bool step(G4Step* aStep, G4TouchableHistory* history)
override;
64 const G4double distance,
67 const G4double stepLength,
68 const G4ThreeVector& mom,
69 const G4ThreeVector& posW,
70 const G4ThreeVector& posIn,
71 const G4ThreeVector& posOut,
72 const G4ThreeVector& posTrack,
77 const G4double hitWeight);
84 G4double
ClosestApproach(G4ThreeVector bwp, G4ThreeVector fwp, G4ThreeVector posIn, G4ThreeVector posOut,
85 G4ThreeVector& hitPosition, G4ThreeVector& wirePosition);
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);
128 void GCUBS(
const G4double x,
const G4double y,
const G4double d1,
const G4double d2, G4double a[4]);
145 void Rotat(G4double& x, G4double& y,
152 void Rotat(G4double x[3],
const int mode);
180 void HELWIR(
const G4double xwb4,
192 const G4double B_kG[3],
193 const G4double charge,
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);
241 unsigned short areNeighbors(
unsigned short iCLayer,
unsigned short iSuperLayer,
unsigned short iLayer,
unsigned short iWire,
242 const WireID& otherWireId)
const;
279 std::multimap<unsigned short, CDCSimHit*>
290 const signed short CCW = 1;
291 const signed short CW = -1;
Class to identify a wire inside the CDC.
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.
const signed short CW_IN_NEIGHBOR
Constant for clockwise inwards.
std::vector< int > WireId_in_hit_order(int id0, int id1, int nWires)
Sort wire id.
const signed short CW
Constant for clockwise orientation.
G4double m_thresholdKineticEnergy
Threshold kinetic energy to be stored.
G4int m_magneticField
Magnetic field is on or off.
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.
const signed short CCW
Constant for counterclockwise orientation.
void GCUBS(const G4double x, const G4double y, const G4double d1, const G4double d2, G4double a[4])
void reAssignLeftRightInfo()
Re-assign left/right info.
CDCSensitiveDetector(G4String name, G4double thresholdEnergyDeposit, G4double thresholdKineticEnergy)
Constructor.
void Rotat(G4double &x, G4double &y, G4double &z, const int mode)
Translation method.
bool step(G4Step *aStep, G4TouchableHistory *history) override
Process each step and calculate variables defined in CDCB4VHit.
const signed short CW_OUT_NEIGHBOR
Constant for clockwise outwards.
std::vector< CDCSimHit * > m_hitWithNegWeight
Vector containing hits with negative weight.
The Class for CDC Geometry Parameters.
Abstract base class for different kinds of events.
int m_hitNumber
The current number of created hits in an event.
CDCGeometryPar * m_cdcgp
Pointer to CDCGeometryPar object.
void for_Rotat(const G4double bfld[3])
Calculates a rotation matrix.
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.
const signed short CCW_OUT_NEIGHBOR
Constant for counterclockwise outwards.
G4int m_nonUniformField
Magnetic field is uniform or non-uniform.
G4double m_thresholdEnergyDeposit
Threshold energy deposit to be stored.
G4double m_minTrackLength
Min.
const signed short CCW_NEIGHBOR
Constant for counterclockwise.
void setModifiedLeftRightFlag()
set left/right flag modified for tracking
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.
const signed short CCW_IN_NEIGHBOR
Constant for counterclockwise inwards.
void Initialize(G4HCofThisEvent *) override
Register CDC hits collection into G4HCofThisEvent.
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.
void EndOfEvent(G4HCofThisEvent *) override
Do what you want to do at the beginning of each event (why this is not called ?)
std::multimap< unsigned short, CDCSimHit * > m_hitWithPosWeight
Map containing hits with positive weight.
~CDCSensitiveDetector()
Destructor.
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).
G4bool m_wireSag
Switch to activate wire sag effect.
const signed short CW_NEIGHBOR
Constant for clockwise.
G4bool m_modifiedLeftRightFlag
Switch for left/right flag modified for tracking.
G4double m_brot[3][3]
a rotation matrix.