ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
doiPETAnalysis.hh
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file doiPETAnalysis.hh
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 
26 //GEANT4 - Depth-of-Interaction enabled Positron emission tomography (PET) advanced example
27 
28 //Authors and contributors
29 
30 // Author list to be updated, with names of co-authors and contributors from National Institute of Radiological Sciences (NIRS)
31 
32 // Abdella M. Ahmed (1, 2), Andrew Chacon (1, 2), Harley Rutherford (1, 2),
33 // Hideaki Tashima (3), Go Akamatsu (3), Akram Mohammadi (3), Eiji Yoshida (3), Taiga Yamaya (3)
34 // Susanna Guatelli (2), and Mitra Safavi-Naeini (1, 2)
35 
36 // (1) Australian Nuclear Science and Technology Organisation, Australia
37 // (2) University of Wollongong, Australia
38 // (3) National Institute of Radiological Sciences, Japan
39 
40 
41 #ifndef doiPETAnalysis_h
42 #define doiPETAnalysis_h 1
43 
44 #include "doiPETGlobalParameters.hh"
45 #include "globals.hh"
46 #include <vector>
47 #include <time.h>
48 #include <map>
49 #include <set>
50 #include "G4ThreeVector.hh"
51 #include <iostream>
52 #include <fstream>
53 #include <sstream>
54 #include <iterator>
55 #include <vector>
56 #include <algorithm>
57 
58 #ifdef USEROOT
59 #pragma GCC diagnostic push
60 #pragma GCC diagnostic ignored "-Wshadow"
61 #include "TRandom3.h"
62 #include "TFile.h"
63 #include "TNtuple.h"
64 #include "g4root.hh"
65 #include "TROOT.h"
66 #include "TFile.h"
67 #include "TH1.h"
68 #include "TH2.h"
69 #include "TTree.h"
70 #include "TLeaf.h"
71 #include "TSystem.h"
72 #pragma GCC diagnostic pop
73 #endif
74 
75 class doiPETAnalysisMessenger;
76 
77 //class InteractionInformation;
78 
79 class doiPETAnalysis
80 {
81 private:
82  doiPETAnalysis();
83 
84 public:
85  ~doiPETAnalysis();
86  static doiPETAnalysis* GetInstance();
87  void FindInteractingCrystal();
88  void Open(G4String);
89  void Close();
90  void Delete();
91  void ResetNumberOfHits();
92  void Write(/*G4int, G4int, G4int, G4double*/);
93  void WriteOutput();
94 
95  //void GetIntractionInfomation(InteractionInformation*);
96 
97  void GetParentParticleName(G4String);
98  void GetSizeOfDetector (G4double, G4double, G4double);
99  void GetScatterIndexInPhantom(G4double);
100 
101  void SetSourcePosition(G4ThreeVector);//
102  void SetEventID(G4int);
103 
104  void BlurringParameters();
105  void GetTimeOfAnnihilation(G4double);
106 
107  void PMTPosition();
108  void AngerLogic(G4int, G4int, G4int, G4double, G4double, G4double, G4double);
109  void ReadReflectorPattern();
110 
111  void SetActivity(G4double);
112  void SetIsotopeHalfLife(G4double);
113  void CrystalIDAfterAngerLogic(G4int, G4int, G4int);
114  void TypeOfOutput(G4String);//Single or coincidence list-mode data
115  void CalulateAcquisitionTime();
116  //G4double QuantumEffifciency(G4double);
117  G4double QuantumEffifciency(G4double, G4int, G4int);
118  void ReadOut(G4int, G4int, G4double, G4double, G4ThreeVector, G4double);
119 
120 private:
121  static doiPETAnalysis* instance;
122  doiPETAnalysisMessenger* fAnalysisMessenger;
123  //std::multimap< G4int, InteractionInformation* > mapBlockInteraction;
124  std::set<G4int> setBlockInteraction;
125 
126  G4double upperThreshold, lowerThreshold;
127 
128 
129 
130  //G4ThreeVector sourcePosition;
131 
132 
133 
134  //
135  G4int scatterIndex;
136 
137  G4String parentParticleName;//
138 
139  //
140  G4int numberofInteractions;
141  G4int countCoincidence;
142 
143  G4int numberOfBlocks_total;
144 
145  G4double sizeOfDetector_DOI,sizeOfDetector_axial,sizeOfDetector_tangential;
146 
147  //Virtual position of the PMT
148  G4double signalPMT1, signalPMT2, signalPMT3, signalPMT4;
149 
150  G4double posPMT1x, posPMT2x, posPMT3x, posPMT4x;
151  G4double posPMT1y, posPMT2y, posPMT3y, posPMT4y;
152  G4double posPMT1z, posPMT2z, posPMT3z, posPMT4z;
153 
154  //
155  G4double signalPMT1z, signalPMT2z, signalPMT3z, signalPMT4z;
156  G4double signalPMT1y, signalPMT2y, signalPMT3y, signalPMT4y;
157 
158  //
159  G4double signalZplus, signalZminus;
160  G4double signalYplus, signalYminus;
161  //
162 
163  G4double dist1z, dist2z, dist3z, dist4z, distz;
164  G4double dist1y, dist2y, dist3y, dist4y, disty;
165 
166  G4double shiftCoeff;
167 
168  G4double PositionAngerZ, PositionAngerY;
169 
170  //reflector pattern
171  std::vector<G4int> ireflectorLayer1_Tangential;
172  std::vector<G4int> ireflectorLayer1_Axial;
173  std::vector<G4int> ireflectorLayer2_Tangential;
174  std::vector<G4int> ireflectorLayer2_Axial;
175  std::vector<G4int> ireflectorLayer3_Tangential;
176  std::vector<G4int> ireflectorLayer3_Axial;
177  std::vector<G4int> ireflectorLayer4_Tangential;
178  std::vector<G4int> ireflectorLayer4_Axial;
179  std::vector<G4int> doi_table;
180  //
181 
182  //The number of pixes for the 2D position histogram after Anger Logic calculation
183  G4int numberOfPixel_axial;
184  G4int numberOfPixel_tan;
185 
186  //source position
187  G4double spositionX;
188  G4double spositionY;
189  G4double spositionZ;
190 
191  //interaction position with respect to the crystal axis
192  G4ThreeVector interactionPos;
193 
194 
195  G4double interactionTime;
196 
197  G4int crystalID;//contineous crystal ID in 3D
198  G4int crystalID_2D;
199 
200  G4int prev_eventID;
201 
202  //Single output
203  G4int eventID;
204  G4int blockID;
205  G4int crystalID_axial;
206  G4int crystalID_tangential;
207  G4int DOI_ID;
208  G4double timeStamp;
209  G4double totalEdep;
210 
211  //coincidence output
212  G4int eventID0, eventID1;
213  G4int blockID0, blockID1;
214  G4int crystalID_axial0, crystalID_axial1;
215  G4int crystalID_tangential0, crystalID_tangential1;
216  G4int DOI_ID0, DOI_ID1;
217  G4double timeStamp0, timeStamp1;
218  G4double totalEdep0, totalEdep1;
219 
220  //choice for the user
221  G4bool getSinglesData;
222  G4bool getCoincidenceData;
223 
224  G4String outputData;
225  G4int numberOfHit;
226  std::vector<G4int> eventID_coin;
227  std::vector<G4double> edep_coin;
228  std::vector<G4int>blockID_coin;
229  std::vector<G4int> cryID_axial_coin;
230  std::vector<G4int> cryID_tan_coin;
231  std::vector<G4int> cryDOI_coin;
232  std::vector<G4double> time_coin;
233 
234  //Crystal IDs after Anger Logic calculation
235  G4int crystalIDNew_DOI, crystalIDNew_tan, crystalIDNew_axial;
236 
237  //Crystal ID in the 2D position histogram along the axial and tangetial direction
238  G4int crystalID_in2D_posHist_axial, crystalID_in2D_posHist_tan;
239 
240  //continous crystal ID after after Anger Logic.
241  G4int crystalID_in2D_posHist;
242 
243 
244  //Crystal blurring
245  G4double crystalResolution;
246  G4double crystalResolutionMin;//
247  G4double crystalResolutionMax;//
248 
249  //G4bool variableResolution;
250  G4bool fixedResolution;
251 
252  G4double energyResolution_fixed;
253  std::vector<std::vector<G4double>> energyResolution_cryDependent;
254 
255  G4double crystalEnergyRef;//This 511 keV
256  G4double crystalQuantumEfficiency;//
257  G4double edep_AfterCrystalBlurring;
258  G4double crystalCoeff;
259  G4double sigma_energyResolution;
260 
261  G4double totalTime;
262  G4double prev_totalTime;
263  G4double timeInterval;
264  G4double time_annihil;
265  G4double time_tof;
266  G4double block_DeadTime;
267  G4double module_DeadTime;
268 
269  G4double *blockTime;
270  G4double *moduleTime;
271 
272  //
273  G4double activityNow;
274  G4double InitialActivity;
275  G4double halfLife;
276 
277  G4String simulationType;
278 
279  //
280  //for output file to write results
281  std::ofstream ofs;
282  G4String asciiFileName;
283  G4String rootFileName;
284 
285 #ifdef USEROOT
286  TTree* tSingles;
287  TTree* tCoincidence;
288  //TH1F*hb;
289 #endif
290 
291  //input file to read reflector pattern
292  std::ifstream ifs;
293 };
294 
295 #endif