CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
FiberSD.cc
Go to the documentation of this file.
4 
5 #include "G4VPhysicalVolume.hh"
6 #include "G4PVPlacement.hh"
7 #include "G4HCofThisEvent.hh"
8 #include "G4TouchableHistory.hh"
9 #include "G4Track.hh"
10 #include "G4Step.hh"
11 #include "G4VSolid.hh"
12 #include "G4DynamicParticle.hh"
13 #include "G4ParticleDefinition.hh"
14 #include "G4SDManager.hh"
15 #include "G4ios.hh"
16 
19  const SimTrackManager* manager) :
20  SensitiveCaloDetector(name, cpv, clg, p), theName(name),
21  m_trackManager(manager), theHCID(-1), theHC(0) {
22 
23  collectionName.insert(name);
24  LogDebug("FiberSim") << "***************************************************"
25  << "\n"
26  << "* *"
27  << "\n"
28  << "* Constructing a FiberSD with name " << GetName()
29  << "\n"
30  << "* *"
31  << "\n"
32  << "***************************************************";
33  theShower = new HFShower(name, cpv, p, 1);
34 
35  //
36  // Now attach the right detectors (LogicalVolumes) to me
37  //
38  std::vector<std::string> lvNames = clg.logicalNames(name);
39  this->Register();
40  for (std::vector<std::string>::iterator it=lvNames.begin();
41  it !=lvNames.end(); it++){
42  this->AssignSD(*it);
43  LogDebug("FiberSim") << "FiberSD : Assigns SD to LV " << (*it);
44  }
45 }
46 
48 
49  if (theShower) delete theShower;
50  if (theHC) delete theHC;
51 }
52 
53 void FiberSD::Initialize(G4HCofThisEvent * HCE) {
54 
55  LogDebug("FiberSim") << "FiberSD : Initialize called for " << GetName();
56  theHC = new FiberG4HitsCollection(GetName(), collectionName[0]);
57  if (theHCID<0)
58  theHCID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
59  HCE->AddHitsCollection(theHCID, theHC);
60 
61 }
62 
63 G4bool FiberSD::ProcessHits(G4Step * aStep, G4TouchableHistory*) {
64 
65  //std::vector<HFShower::Hit> hits = theShower->getHits(aStep);
66  double zoffset = 1000;
67  std::vector<HFShower::Hit> hits = theShower->getHits(aStep,true,zoffset);
68 
69 
70  if (hits.size() > 0) {
71  std::vector<HFShowerPhoton> thePE;
72  for (unsigned int i=0; i<hits.size(); i++) {
73  //std::cout<<"hit position z "<<hits[i].position.z()<<std::endl;
74  HFShowerPhoton pe = HFShowerPhoton(hits[i].position.x(),
75  hits[i].position.y(),
76  hits[i].position.z(),
77  hits[i].wavelength, hits[i].time);
78  thePE.push_back(pe);
79  }
80  int trackID = aStep->GetTrack()->GetTrackID();
81  G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
82  const G4VTouchable* touch = preStepPoint->GetTouchable();
83  G4LogicalVolume* lv = touch->GetVolume(0)->GetLogicalVolume();
84  int depth = (touch->GetReplicaNumber(0))%10;
85  int detID = setDetUnitId(aStep);
86  math::XYZPoint theHitPos(preStepPoint->GetPosition().x(),
87  preStepPoint->GetPosition().y(),
88  preStepPoint->GetPosition().z());
89  //std::cout<<"presteppoint position z "<<preStepPoint->GetPosition().z()<<std::endl;
90 
91  FiberG4Hit *aHit = new FiberG4Hit(lv, detID, depth, trackID);
92  std::cout<<"hit size "<<hits.size()<<" npe"<<aHit->npe()<<std::endl;
93  std::cout<<"pre hit position "<<aHit->hitPos()<<std::endl;
94  aHit->setNpe(hits.size());
95  aHit->setPos(theHitPos);
96  aHit->setTime(preStepPoint->GetGlobalTime());
97  aHit->setPhoton(thePE);
98  std::cout<<"ShowerPhoton position "<<thePE[0].x()<<" "<<thePE[0].y()<<" "<<thePE[0].z()<<std::endl;
99 
100  LogDebug("FiberSim") << "FiberSD: Hit created at " << lv->GetName()
101  << " DetID: " << aHit->towerId() << " Depth: "
102  << aHit->depth() << " Track ID: " << aHit->trackId()
103  << " Nb. of Cerenkov Photons: " << aHit->npe()
104  << " Time: " << aHit->time() << " at "
105  << aHit->hitPos();
106  for (unsigned int i=0; i<thePE.size(); i++)
107  LogDebug("FiberSim") << "FiberSD: PE[" << i << "] " << thePE[i];
108 
109  theHC->insert(aHit);
110  }
111  return true;
112 }
113 
114 void FiberSD::EndOfEvent(G4HCofThisEvent * HCE) {
115 
116  LogDebug("FiberSim") << "FiberSD: Sees" << theHC->entries() << " hits";
117  clear();
118  std::cout<<"theHC entries = "<<theHC->entries()<<std::endl;
119 }
120 
121 void FiberSD::clear() {}
122 
124 
126 
127 void FiberSD::update(const BeginOfRun *) {}
128 
130 
131 void FiberSD::update(const ::EndOfEvent *) {}
132 
134 
135 uint32_t FiberSD::setDetUnitId(G4Step* aStep) {
136  const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
137  int fibre = (touch->GetReplicaNumber(1))%10;
138  int cell = (touch->GetReplicaNumber(2));
139  int tower = (touch->GetReplicaNumber(3));
140  return ((tower*1000+cell)*10+fibre);
141 }
142 
#define LogDebug(id)
int i
Definition: DBlmapReader.cc:9
std::vector< PCaloHit > PCaloHitContainer
std::vector< std::string > logicalNames(std::string &readoutName)
virtual void DrawAll()
Definition: FiberSD.cc:123
type of data representation of DDCompactView
Definition: DDCompactView.h:77
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
virtual void update(const BeginOfRun *)
This routine will be called when the appropriate signal arrives.
Definition: FiberSD.cc:127
virtual void fillHits(edm::PCaloHitContainer &, std::string)
Definition: FiberSD.cc:143
std::string const collectionName[nCollections]
Definition: Collections.h:45
virtual void PrintAll()
Definition: FiberSD.cc:125
std::vector< Hit > getHits(G4Step *aStep, double weight)
Definition: HFShower.cc:71
FiberSD(std::string, const DDCompactView &, SensitiveDetectorCatalog &, edm::ParameterSet const &, const SimTrackManager *)
Definition: FiberSD.cc:17
HFShower * theShower
Definition: FiberSD.h:57
virtual uint32_t setDetUnitId(G4Step *)
Definition: FiberSD.cc:135
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
FiberG4HitsCollection * theHC
Definition: FiberSD.h:60
G4THitsCollection< FiberG4Hit > FiberG4HitsCollection
Definition: FiberG4Hit.h:59
virtual void AssignSD(std::string &vname)
virtual void clearHits()
Definition: FiberSD.cc:133
tuple cout
Definition: gather_cfg.py:121
virtual void Initialize(G4HCofThisEvent *HCE)
Definition: FiberSD.cc:53
virtual void EndOfEvent(G4HCofThisEvent *HCE)
Definition: FiberSD.cc:114
virtual ~FiberSD()
Definition: FiberSD.cc:47
virtual void clear()
Definition: FiberSD.cc:121
G4int theHCID
Definition: FiberSD.h:59
virtual G4bool ProcessHits(G4Step *aStep, G4TouchableHistory *ROhist)
Definition: FiberSD.cc:63