1515 Includes functionalities from:
1616 GateSpblurring
1717 GateCC3DlocalSpblurring
18- GateDoIModels (TODO)
18+ GateDoIModels
1919
2020 Previous authors: Steven.Staelens@rug.ac.be(?), AE
2121
@@ -65,12 +65,12 @@ GateSpatialResolution::GateSpatialResolution(GateSinglesDigitizer *digitizer, G4
6565 m_fwhmXDistrib(0 ),
6666 m_fwhmYDistrib(0 ),
6767 m_fwhmZDistrib(0 ),
68- m_nameAxis(" YZ " ),
68+ m_nameAxis(" " ),
6969 m_fwhmXDistrib2D(0 ),
7070 m_fwhmYDistrib2D(0 ),
7171 m_fwhmZDistrib2D(0 ),
7272 m_IsConfined(true ),
73- m_UseTruncatedGaussian(true ),
73+ m_UseTruncatedGaussian(false ),
7474 m_Navigator(0 ),
7575 m_Touchable(0 ),
7676 m_systemDepth(-1 ),
@@ -129,12 +129,30 @@ void GateSpatialResolution::SetSpatialResolutionParameters() {
129129 abort ();
130130 }
131131
132+
133+ if (m_nameAxis == " YZ" && (m_fwhmYDistrib2D ==0 || m_fwhmZDistrib2D ==0 ))
134+ {
135+ G4cout << " ***ERROR*** Spatial Resolution 2D distribution: you choose YZ axis but one or both parameters /fwhmYDistrib2D or /m_fwhmZDistrib2D are not defined" << G4endl;
136+ abort ();
137+ }
138+ if (m_nameAxis == " XZ" && (m_fwhmXDistrib2D ==0 || m_fwhmZDistrib2D ==0 ))
139+ {
140+ G4cout << " ***ERROR*** Spatial Resolution 2D distribution: you choose XZ axis but one or both parameters /fwhmXDistrib2D or /m_fwhmZDistrib2D are not defined" << G4endl;
141+ abort ();
142+ }
143+
144+
145+
146+
147+
148+
132149}
133150
134151
135152
136153void GateSpatialResolution::Digitize (){
137154
155+
138156 GateVSystem* m_system = ((GateSinglesDigitizer*)this ->GetDigitizer ())->GetSystem ();
139157
140158 if (m_IsFirstEntrance) {
@@ -190,6 +208,17 @@ void GateSpatialResolution::Digitize(){
190208 GateDigi* inputDigi;
191209
192210
211+ if (!m_IsConfined && !m_Navigator)
212+ {
213+ // Getting world Volume
214+ // Do not use from TransportationManager as it is not recommended
215+ G4Navigator *navigator = G4TransportationManager::GetTransportationManager ()->GetNavigatorForTracking ();
216+ G4VPhysicalVolume *WorldVolume = navigator->GetWorldVolume ();
217+ m_Navigator = new G4Navigator ();
218+ m_Navigator->SetWorldVolume (WorldVolume);
219+ }
220+
221+
193222 if (IDC)
194223 {
195224 G4int n_digi = IDC->entries ();
@@ -208,39 +237,41 @@ void GateSpatialResolution::Digitize(){
208237 G4double Pz = P.z ();
209238 G4double stddevX = 0 ., stddevY = 0 ., stddevZ = 0 .;
210239
211- // Use the configured axis pair (m_nameAxis) to evaluate Value2D.
212240 // Allowed pairs for PET: "XZ" or "YZ" (default "YZ").
213- if (m_fwhmXDistrib2D || m_fwhmYDistrib2D || m_fwhmZDistrib2D) {
214- if (m_nameAxis == " XZ" ) {
215- if (m_fwhmXDistrib2D) stddevX = m_fwhmXDistrib2D->Value2D (P.x () * mm, P.z () * mm);
216- if (m_fwhmYDistrib2D) stddevY = m_fwhmYDistrib2D->Value2D (P.x () * mm, P.z () * mm);
217- if (m_fwhmZDistrib2D) stddevZ = m_fwhmZDistrib2D->Value2D (P.x () * mm, P.z () * mm);
218- } else { // YZ
219- if (m_fwhmXDistrib2D) stddevX = m_fwhmXDistrib2D->Value2D (P.y () * mm, P.z () * mm);
220- if (m_fwhmYDistrib2D) stddevY = m_fwhmYDistrib2D->Value2D (P.y () * mm, P.z () * mm);
221- if (m_fwhmZDistrib2D) stddevZ = m_fwhmZDistrib2D->Value2D (P.y () * mm, P.z () * mm);
222- }
241+ if (m_nameAxis == " XZ" ) { // if 2D distributions are defined
242+ stddevX = m_fwhmXDistrib2D->Value2D (P.x () * mm, P.z () * mm);
243+ stddevZ = m_fwhmZDistrib2D->Value2D (P.x () * mm, P.z () * mm);
223244 }
224- else {
225-
226- if (m_fwhmXDistrib) stddevX = m_fwhmXDistrib->Value (P.x () * mm);
227- else if (fwhmX) stddevX = fwhmX / GateConstants::fwhm_to_sigma;
228-
229- if (m_fwhmYDistrib) stddevY = m_fwhmYDistrib->Value (P.y () * mm);
230- else if (fwhmY) stddevY = fwhmY / GateConstants::fwhm_to_sigma;
231-
232- if (m_fwhmZDistrib) stddevZ = m_fwhmZDistrib->Value (P.z () * mm);
233- else if (fwhmZ) stddevZ = fwhmZ / GateConstants::fwhm_to_sigma;
234-
245+ if (m_nameAxis == " YZ" ) { // YZ
246+ stddevY = m_fwhmYDistrib2D->Value2D (P.y () * mm, P.z () * mm);
247+ stddevZ = m_fwhmZDistrib2D->Value2D (P.y () * mm, P.z () * mm);
235248 }
236-
237-
238-
249+
250+
251+ if (m_fwhmXDistrib) stddevX = m_fwhmXDistrib->Value (P.x () * mm);
252+ else if (fwhmX)
253+ stddevX = fwhmX / GateConstants::fwhm_to_sigma;
254+
255+
256+ if (m_fwhmYDistrib) stddevY = m_fwhmYDistrib->Value (P.y () * mm);
257+ else if (fwhmY) stddevY = fwhmY / GateConstants::fwhm_to_sigma;
258+
259+ if (m_fwhmZDistrib) stddevZ = m_fwhmZDistrib->Value (P.z () * mm);
260+ else if (fwhmZ) stddevZ = fwhmZ / GateConstants::fwhm_to_sigma;
261+
262+ // }
263+
264+
265+ if (!stddevX || !stddevY || !stddevZ)
266+ {
267+ G4cout << " *** ERROR*** GateSpatialResolution::Digitize. Resolution for one of the axis is not set. If you really don't want to apply spatial resolution on one of the axis, please, contact kochebina [a] gmail (dot) com " << G4endl;
268+ abort ();
269+ }
239270 G4double PxNew ;
240271 G4double PyNew ;
241272 G4double PzNew ;
242273
243- // store the computed stddevs into the digi for later ROOT output
274+ // store the computed stddevs into the digi for later ROOT output
244275 m_outputDigi->SetSpatialRes2DStdDevX (stddevX);
245276 m_outputDigi->SetSpatialRes2DStdDevY (stddevY);
246277 m_outputDigi->SetSpatialRes2DStdDevZ (stddevZ);
@@ -314,15 +345,11 @@ void GateSpatialResolution::Digitize(){
314345 if (PxNew>Xmax) PxNew=Xmax;
315346 if (PyNew>Ymax) PyNew=Ymax;
316347 if (PzNew>Zmax) PzNew=Zmax;
348+
317349 m_outputDigi->SetLocalPos (G4ThreeVector (PxNew,PyNew,PzNew)); // TC
318- m_outputDigi->SetGlobalPos (m_outputDigi->GetVolumeID ().MoveToAncestorVolumeFrame (m_outputDigi->GetLocalPos ())); // TC
319350
320- // Getting world Volume
321- // Do not use from TransportationManager as it is not recommended
322- G4Navigator *navigator = G4TransportationManager::GetTransportationManager ()->GetNavigatorForTracking ();
323- G4VPhysicalVolume *WorldVolume = navigator->GetWorldVolume ();
324- m_Navigator = new G4Navigator ();
325- m_Navigator->SetWorldVolume (WorldVolume);
351+
352+ m_outputDigi->SetGlobalPos (m_outputDigi->GetVolumeID ().MoveToAncestorVolumeFrame (m_outputDigi->GetLocalPos ())); // TC
326353
327354 G4VPhysicalVolume* PV = m_Navigator->LocateGlobalPointAndSetup (m_outputDigi->GetGlobalPos ());
328355 m_Touchable = m_Navigator->CreateTouchableHistoryHandle ();
@@ -332,6 +359,7 @@ void GateSpatialResolution::Digitize(){
332359 {
333360 UpdateVolumeID ();
334361 }
362+
335363 m_OutputDigiCollection->insert (m_outputDigi);
336364
337365
0 commit comments