doiPET/src/doiPETSteppingAction.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025
0026 //GEANT4 - Depth-of-Interaction enabled Positron emission tomography (PET) advanced example
0027
0028 //Authors and contributors
0029
0030 // Author list to be updated, with names of co-authors and contributors from National Institute of Radiological Sciences (NIRS)
0031
0032 // Abdella M. Ahmed (1, 2), Andrew Chacon (1, 2), Harley Rutherford (1, 2),
0033 // Hideaki Tashima (3), Go Akamatsu (3), Akram Mohammadi (3), Eiji Yoshida (3), Taiga Yamaya (3)
0034 // Susanna Guatelli (2), and Mitra Safavi-Naeini (1, 2)
0035
0036 // (1) Australian Nuclear Science and Technology Organisation, Australia
0037 // (2) University of Wollongong, Australia
0038 // (3) National Institute of Radiological Sciences, Japan
0039
0040
0041 #include "doiPETSteppingAction.hh"
0042 #include "doiPETAnalysis.hh"
0043 #include "doiPETRun.hh"
0044 #include "doiPETEventAction.hh"
0045
0046 #include "G4Step.hh"
0047 #include "G4RunManager.hh"
0048 #include "G4SteppingManager.hh"
0049
0050 #include "G4SystemOfUnits.hh"
0051 #include "G4PhysicalConstants.hh"
0052 #include <limits>
0053
0054 ////////// Constructor //////////////////////////////////////////////
0055 /*doiPETSteppingAction::doiPETSteppingAction()
0056 {;}*/
0057 doiPETSteppingAction::doiPETSteppingAction()
0058     : G4UserSteppingAction()
0059 { }
0060
0061 ///////// Destructor ////////////////////////////////////////////////
0062 doiPETSteppingAction::~doiPETSteppingAction()
0063 {;}
0064
0065 ///////// UserSteppingAction ///////////////////////////////////////
0066 void doiPETSteppingAction::UserSteppingAction(const G4Step* aStep)
0067 {
0068     doiPETRun* run = static_cast<doiPETRun*>(G4RunManager::GetRunManager()->GetNonConstCurrentRun());
0069
0070     G4Track* track = aStep->GetTrack();
0071     G4String volumeName = aStep->GetPreStepPoint()->GetPhysicalVolume()->GetName();
0072     G4double edep = aStep->GetTotalEnergyDeposit();
0073     G4String processName =aStep->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName();
0074     G4String particleName = track->GetDynamicParticle()->GetDefinition()->GetParticleName();
0075     //G4int eventID = G4RunManager::GetRunManager()->GetCurrentEvent() -> GetEventID();
0076     G4ThreeVector pos = track->GetPosition();
0077
0078
0079     //G4StepPoint* p2 = aStep->GetPostStepPoint();
0080     G4StepPoint* p1 = aStep->GetPreStepPoint();
0081     G4ThreeVector coord1 = p1->GetPosition();
0082
0083     //The following is to get the local position of the interaction with respect to the crystal block volume
0084     const G4AffineTransform transformation =  p1->GetTouchable()-> GetHistory()->GetTransform(2);
0085     G4ThreeVector localPosition = transformation.TransformPoint(coord1);
0086
0087     G4int blockID;
0088     G4int crystalID;
0089
0090     G4int scatterIndex = 0;//For checking
0091
0092     //If an event is created in a cold region or outside the physical volume, then it will be killed.
0093     if((volumeName != "phantom_physicalV" || volumeName == "coldRegion_physicalV") && particleName== "e+"){
0094         track->SetTrackStatus(fStopAndKill);
0095         return;
0096     }
0097
0098     //Get the source position if the process name is annihilation
0099     if(processName == "annihil" && volumeName=="phantom_physicalV"){
0100         doiPETAnalysis::GetInstance()->SetSourcePosition(pos);
0101         run->SetAnnihilationTime((track->GetGlobalTime()));
0102     }
0103
0104     //get event ID
0105     //doiPETAnalysis::GetInstance()->SetEventID (eventID);
0106
0107     //Get scatter information in the phantom by the annihilation photon before detected by the detector. Note that the scatter index is initialized to 0.
0108     //If there energy deposition in the phantom by the photon, then there is scatter, the index is 1, and if not it is 0.
0109     if(edep > 0 && particleName == "gamma"  && volumeName == "phantom_physicalV"){
0110         scatterIndex = 1;
0111         doiPETAnalysis::GetInstance()->SetScatterIndexInPhantom(scatterIndex);
0112     }
0113
0114     ///////////////////    Retrive (Extract) information in the crystal ///////////////////////////
0115     if(edep>0. && volumeName=="Crystal_physicalV"){
0116
0117         //get the copy number of the block
0118         blockID = aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(2);
0119
0120         //get the crystal copy number
0121         crystalID = aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(0);
0122
0123         //G4cout<<localPosition.x()<<" "<<localPosition.y()<<" "<<localPosition.z()<<G4endl;
0124
0125         //get event ID
0126         //doiPETAnalysis::GetInstance()->GetEventID (eventID);
0127
0128         //define a pointer for extracting information
0129         InteractionInformation* ExtractIntInfo = new InteractionInformation();
0130
0131         //get the energy deposition of each interaction
0132         ExtractIntInfo->SetEdep( edep );
0133
0134
0135         //get the blockID
0136         ExtractIntInfo->SetBlockNo( blockID );
0137
0138
0139         //get the crystal ID
0140         ExtractIntInfo->SetCrystalNo( crystalID );
0141
0142
0143         //get local position interaction in the crystal
0144         ExtractIntInfo->SetInteractionPositionInCrystal(localPosition);
0145
0146         //get the global time of the interaction with the crystal
0147         ExtractIntInfo->SetGlobalTime( track->GetGlobalTime());
0148
0149         //pass all the obtained information to the doiPETAnalysis class
0150         //doiPETAnalysis::GetInstance()->GetIntractionInfomation(ExtractIntInfo);
0151         run->GetIntractionInfomation(ExtractIntInfo);
0152
0153     }
0154 }