QwGeant4
QweakSimVDC_DriftCellFrontSD.cc
Go to the documentation of this file.
1 //=============================================================================
2 //
3 // ---------------------------
4 // | Doxygen File Information |
5 // ---------------------------
6 //
7 /**
8 
9  \file QweakSimVDC_DriftCellFrontSD.cc
10 
11  $Revision: 1.5 $
12  $Date: 2006/05/05 21:47:42 $
13 
14  \author Klaus Hans Grimm
15 
16 */
17 //=============================================================================
18 
19 //=============================================================================
20 // -----------------------
21 // | CVS File Information |
22 // -----------------------
23 //
24 // Last Update: $Author: grimm $
25 // Update Date: $Date: 2006/05/05 21:47:42 $
26 // CVS/RCS Revision: $Revision: 1.5 $
27 // Status: $State: Exp $
28 //
29 // ===================================
30 // CVS Revision Log at end of file !!
31 // ===================================
32 //
33 //============================================================================
34 
35 
36 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
37 
39 
40 // user includes
44 
45 
46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
47 
53 
54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56 :G4VSensitiveDetector(name)
57 {
58  //G4cout << G4endl << "###### Calling QweakSimVDC_DriftCellFrontSD::QweakSimVDC_DriftCellFrontSD() " << G4endl << G4endl;
59 
60  collectionName.insert("DriftCellFrontCollection");
61 
62  DC_ID = -1;
63 
64  //G4cout << G4endl << "###### Leaving QweakSimVDC_DriftCellFrontSD::QweakSimVDC_DriftCellFrontSD() " << G4endl << G4endl;
65 }
66 
67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
69 {
70  //G4cout << G4endl << "###### Calling QweakSimVDC_DriftCellFrontSD::~QweakSimVDC_DriftCellFrontSD() " << G4endl << G4endl;
71 
72 }
73 
74 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
75 void QweakSimVDC_DriftCellFrontSD::Initialize(G4HCofThisEvent* HCE)
76 {
77  //G4cout << G4endl << "###### Calling QweakSimVDC_DriftCellFrontSD::Initialize() " << G4endl << G4endl;
78 
79  DC_hitsCollection = new QweakSimVDC_DriftCellHitsCollection(SensitiveDetectorName,collectionName[0]);
80 
81  // for(int i=0;i<DCNumberPerPlane;i++)
82  //{ DC_hitsCollection->insert(new QweakSimVDC_DriftCellHit); }
83 
84  if(DC_ID<0)
85  { DC_ID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]); }
86  //{ DC_ID = G4SDManager::GetSDMpointer()->GetCollectionID(DC_hitsCollection); }
87  HCE->AddHitsCollection( DC_ID, DC_hitsCollection );
88 
89  //G4cout << G4endl << "###### Leaving QweakSimVDC_DriftCellFrontSD::Initialize() " << G4endl << G4endl;
90 }
91 
92 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
93 G4bool QweakSimVDC_DriftCellFrontSD::ProcessHits(G4Step* aStep,G4TouchableHistory* )
94 {
95  //G4cout << G4endl << "###### Calling QweakSimVDC_DriftCellFrontSD::ProcessHits() " << G4endl << G4endl;
96 
97  //G4double charge = aStep->GetTrack()->GetDefinition()->GetPDGCharge();
98 
99  // G4cout << "%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Charge =" << charge << G4endl;
100  //
101  // dismiss photons
102  // if(charge==0.) {
103 
104  // G4cout << "==================================================================" << G4endl;
105  // G4cout << ">>>> Neutral Particle detected (e.g. Photon) <<<<<<" << G4endl;
106  // G4cout << ">>>> Up to now we don't count HITs for neutrals, therefore: <<<<<<" << G4endl;
107  // G4cout << ">>>> Aborting QweakSimVDC_DriftCellFrontSD::ProcessHits <<<<<<" << G4endl;
108  // G4cout << "==================================================================" << G4endl;
109 
110  // G4cout << G4endl << "###### Leaving QweakSimVDC_DriftCellFrontSD::ProcessHits() " << G4endl << G4endl;
111  // return false;
112  // }
113 
114 
115 
116 
117  //=====================================================================================
118  G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
119  //G4StepPoint* postStepPoint = aStep->GetPreStepPoint();
120  G4TouchableHistory* theTouchable = (G4TouchableHistory*)(preStepPoint->GetTouchable());
121  G4VPhysicalVolume* physVol = theTouchable->GetVolume();
122  //====================================================================================
123 
124 
125 // G4cout << "==================================================================" << G4endl;
126 // G4cout << ">>>>>> We are in volume :" << physVol->GetName() << "<<<<<<" << G4endl;
127 // G4cout << "==================================================================" << G4endl;
128 
129 
130  if( physVol->GetName().compare("VDC_DriftCellFront_Physical")==0 )
131  {
132  // G4cout << "==================================================================" << G4endl;
133  // G4cout << ">>>>>>> Particle crossing : VDC_DriftCellFront_Physical <<<<<<<<<" << G4endl;
134  // G4cout << "==================================================================" << G4endl;
135 
136  //const G4VProcess* process = aStep->GetPostStepPoint()->GetProcessDefinedStep();
137  //G4cout << "@@@@@@@@@@@@@> Physical Processes in VDC_DriftCellFront_Physical = "<< process->GetProcessName() << G4endl;
138 
139  }
140  else
141  {
142  // G4cout << "==================================================================" << G4endl;
143  // G4cout << ">>>> ERROR: We are not within a physical DriftCell Volume <<<<<<" << G4endl;
144  // G4cout << ">>>> Aborting QweakSimVDC_DriftCellFrontSD::ProcessHits <<<<<<" << G4endl;
145  // G4cout << "==================================================================" << G4endl;
146 
147  return false;
148  }
149 
150  //================================================================================================
151  // here we are only interested if the particle is crossing a boundary
152  // otherwise we get additional hits due to physical processes within the sensitive volume
153  //
154  // If you want to identify the first step in a volume: pick fGeomBoundary status in preStepPoint
155  // If you want to identify a step out of a volume : pick fGeomBoundary status in postStepPoint
156  //=================================================================================================
157  //
158  if (preStepPoint->GetStepStatus() != fGeomBoundary)
159  {
160  // G4cout << "=======================================================================" << G4endl;
161  // G4cout << ">>>> QweakSimVDC_DriftCellFrontSD: We are NOT crossing a boundary <<<<<<" << G4endl;
162  // G4cout << ">>>> Aborting QweakSimVDC_DriftCellFrontSD::ProcessHits <<<<<<" << G4endl;
163  // G4cout << "=======================================================================" << G4endl;
164  //
165  // G4cout << G4endl << "###### Leaving QweakSimVDC_DriftCellFrontSD::ProcessHits() " << G4endl << G4endl;
166 
167  return false;
168  }
169 
170 
171  G4double currentKineticEnergy = aStep->GetTrack()->GetKineticEnergy();
172  G4double currentTotalEnergy = aStep->GetTrack()->GetTotalEnergy();
173  G4ThreeVector currentMomentumDirection = aStep->GetTrack()->GetMomentumDirection();
174 
175  //G4int trackID = aStep->GetTrack()->GetTrackID();
176  //G4cout << "====> Track ID : " << trackID << G4endl;
177 
178  //G4int parentID = aStep->GetTrack()->GetParentID();
179  //G4cout << "====> Parent ID : " << parentID << G4endl;
180 
181  // Vertex: where this track was generated
182  G4ThreeVector originVertexPosition = aStep->GetTrack()->GetVertexPosition();
183  G4ThreeVector originVertexMomentumDirection = aStep->GetTrack()->GetVertexMomentumDirection();
184  G4double originVertexKineticEnergy = aStep->GetTrack()->GetVertexKineticEnergy();
185 
186 // G4cout << "====> originVertexPosition_X [cm] : " << originVertexPosition.x()/cm << G4endl;
187 // G4cout << "====> originVertexPosition_Y [cm] : " << originVertexPosition.y()/cm << G4endl;
188 // G4cout << "====> originVertexPosition_Z [cm] : " << originVertexPosition.z()/cm << G4endl;
189 
190 // G4cout << "====> originVertexMomentumDirection_X : " << originVertexMomentumDirection.x() << G4endl;
191 // G4cout << "====> originVertexMomentumDirection_Y : " << originVertexMomentumDirection.y() << G4endl;
192 // G4cout << "====> originVertexMomentumDirection_Z : " << originVertexMomentumDirection.z() << G4endl;
193 
194 // G4cout << "====> originVertexKineticEnergy [Mev] : " << originVertexKineticEnergy/MeV << G4endl;
195 
196 
197  //G4int MotherCopyNo = theTouchable->GetVolume(1)->GetCopyNo(); // one Mother Volume
198  //G4int DriftCellCopyNo = theTouchable->GetVolume()->GetCopyNo(); // but several Driftcells per MV
199  G4int DriftCellReplicaNo = theTouchable->GetReplicaNumber(); // but several Driftcells per MV
200 
201 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront MV CopyNumber :" << MotherCopyNo << G4endl;
202 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront CopyNumber :" << DriftCellCopyNo << G4endl;
203 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront ReplicaNumber :" << DriftCellReplicaNo << G4endl;
204 
205 
206  G4ThreeVector worldPos = preStepPoint->GetPosition();
207 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront WorldPos.X [cm]:" << worldPos.x()/cm << G4endl;
208 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront WorldPos.Y [cm]:" << worldPos.y()/cm << G4endl;
209 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront WorldPos.Z [cm]:" << worldPos.z()/cm << G4endl;
210 
211  G4ThreeVector localPos = theTouchable->GetHistory()->GetTopTransform().TransformPoint(worldPos);
212 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront LocalPos.X [cm]:" << localPos.x()/cm << G4endl;
213 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront LocalPos.Y [cm]:" << localPos.y()/cm << G4endl;
214 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% VDC_DriftCellFront LocalPos.Z [cm]:" << localPos.z()/cm << G4endl;
215 
216  QweakSimVDC_DriftCellHit* aHit = new QweakSimVDC_DriftCellHit(DriftCellReplicaNo);
217  //QweakSimVDC_DriftCellHit* aHit = (*DC_hitsCollection)[DriftCellReplicaNo];
218 
219 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% got QweakSimVDC_DriftCellHit* aHit" << G4endl;
220 
221  aHit->StoreWorldPos(worldPos);
222  aHit->StoreLocalPos(localPos);
223  aHit->StoreTime(preStepPoint->GetGlobalTime());
224 
225  aHit->StoreOriginVertexPosition(originVertexPosition);
226  aHit->StoreOriginVertexKineticEnergy(originVertexKineticEnergy);
227  aHit->StoreOriginVertexMomentumDirection(originVertexMomentumDirection);
228 
229  aHit->StoreMomentumDirection(currentMomentumDirection);
230  aHit->StoreKineticEnergy(currentKineticEnergy);
231  aHit->StoreTotalEnergy(currentTotalEnergy);
232 
233 
234  G4double myDCUPlaneWireAngle = 90.0*degree - DCUPlaneWireAngle;
235  G4double myDCVPlaneWireAngle = -90.0*degree - DCVPlaneWireAngle;
236 
237 // G4cout << " DCWidthOnFrame [mm] = " << DCWidthOnFrame/mm << G4endl;
238 // G4cout << " DCFullThickness [mm] = " << DCFullThickness/mm << G4endl;
239 // G4cout << " DCUPlaneWireAngle [deg] = " << myDCUPlaneWireAngle/degree << G4endl;
240 // G4cout << " DCVPlaneWireAngle [deg] = " << myDCVPlaneWireAngle/degree << G4endl;
241 
244  aHit->StoreDCUPlaneWireAngle(myDCUPlaneWireAngle);
245  aHit->StoreDCVPlaneWireAngle(myDCVPlaneWireAngle);
246 
247 
248 
249  // check if it is first touch
250  if(!(aHit->GetLogV()))
251  {
252  // store translation and rotation matrix of the drift cell
253  // for the sake of drawing the hit
254  aHit->StoreLogV(physVol->GetLogicalVolume());
255  G4AffineTransform aTrans = theTouchable->GetHistory()->GetTopTransform();
256  aTrans.Invert();
257  aHit->StoreCellRot(aTrans.NetRotation());
258  aHit->StoreCellPos(aTrans.NetTranslation());
259  }
260 
261  if( physVol->GetName().compare("VDC_DriftCellFront_Physical") == 0 ) { aHit->StoreDriftCellPlaneID(0);}
262  if( physVol->GetName().compare("VDC_DriftCellBack_Physical") == 0 ) { aHit->StoreDriftCellPlaneID(1);}
263 
264 // G4cout << "%%%%%%%%%%%%%%%%%%%%%%% before inserting the hit" << G4endl;
265 
266  DC_hitsCollection->insert(aHit);
267 
268 // G4cout << G4endl << "###### Leaving QweakSimVDC_DriftCellFrontSD::ProcessHits() " << G4endl << G4endl;
269 
270  return true;
271 }
272 
273 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
274 
275 //void QweakSimVDC_DriftCellFrontSD::EndOfEvent(G4HCofThisEvent* HCE)
277 {
278  //G4cout << G4endl << "###### Calling QweakSimVDC_DriftCellFrontSD::EndOfEvent() " << G4endl << G4endl;
279 
280 // G4int NbDCHits = DC_hitsCollection->entries();
281 
282 // G4cout << "\n-------->Hits Collection: in this event they are " << NbDCHits
283 // << " hits in the Drift Cells : " << G4endl;
284 // for (G4int i=0;i<NbDCHits;i++) (*DC_hitsCollection)[i]->Print();
285 
286 
287 
288  //G4cout << G4endl << "###### Leaving QweakSimVDC_DriftCellFrontSD::EndOfEvent() " << G4endl << G4endl;
289 
290 }
291 
292 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
293 
294 //=======================================================
295 // -----------------------
296 // | CVS File Information |
297 // -----------------------
298 //
299 // $Revisions$
300 // $Log: QweakSimVDC_DriftCellFrontSD.cc,v $
301 // Revision 1.5 2006/05/05 21:47:42 grimm
302 // Seperation into DriftCellBack and DriftCellFront for testing purpose.
303 //
304 // Revision 1.4 2006/01/06 21:41:29 grimm
305 // Blank line added. No changes
306 //
307 // Revision 1.3 2006/01/02 13:24:31 grimm
308 // Commenting the originVertexPosition cut (z> 635*cm). Now I'm allowing to record all charged hits of the drift cell regardless of the vertex origin.
309 // (With the target cut I got ~20k drift cell hits out of 300k primary events. The wire planes without origin cut get a hit quote of 30%)
310 //
311 // Revision 1.2 2005/12/27 19:19:23 grimm
312 // - Redesign of Doxygen header containing CVS info like revision and date
313 // - Added CVS revision log at the end of file
314 //
315 //
void StoreWorldPos(G4ThreeVector xyz)
void StoreDCVPlaneWireAngle(G4double dca)
void StoreCellRot(G4RotationMatrix rmat)
G4THitsCollection< QweakSimVDC_DriftCellHit > QweakSimVDC_DriftCellHitsCollection
QweakSimVDC_DriftCellHitsCollection * DC_hitsCollection
void StoreCellPos(G4ThreeVector xyz)
void StoreOriginVertexKineticEnergy(G4double ekin)
const G4LogicalVolume * GetLogV() const
void StoreOriginVertexPosition(G4ThreeVector xyz)
void StoreLocalPos(G4ThreeVector xyz)
void StoreOriginVertexMomentumDirection(G4ThreeVector pxyz)
void StoreKineticEnergy(G4double ekin)
Region 3 Vertical Drift Chamber Drift Cell Hit.
void StoreDriftCellPlaneID(G4int dcplane_id)
void StoreLogV(G4LogicalVolume *val)
G4bool ProcessHits(G4Step *aStep, G4TouchableHistory *ROhist)
void StoreMomentumDirection(G4ThreeVector pxyz)
void StoreDCUPlaneWireAngle(G4double dca)