microbeam/src/SteppingAction.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 // This example is provided by the Geant4-DNA collaboration
0027 // Any report or published results obtained using the Geant4-DNA software
0028 // shall cite the following Geant4-DNA collaboration publication:
0029 // Med. Phys. 37 (2010) 4692-4708
0030 // The Geant4-DNA web site is available at http://geant4-dna.org
0031 //
0032 // If you use this example, please cite the following publication:
0033 // Rad. Prot. Dos. 133 (2009) 2-11
0034 //
0035 #include "SteppingAction.hh"
0036 #include "RunAction.hh"
0037 #include "DetectorConstruction.hh"
0038 #include "G4AnalysisManager.hh"
0039 #include "G4SystemOfUnits.hh"
0040 #include "G4SteppingManager.hh"
0041 #include "G4Alpha.hh"
0042 #include "G4Electron.hh"
0043
0044 SteppingAction::SteppingAction(RunAction* run,const DetectorConstruction* det)
0045 :fRun(run),fDetector(det)
0046 {}
0047
0048 SteppingAction::~SteppingAction()
0049 {}
0050
0051 void SteppingAction::UserSteppingAction(const G4Step* aStep)
0052
0053 {
0054   // Analysis manager
0055
0056   G4AnalysisManager* man = G4AnalysisManager::Instance();
0057
0058   // Read phantom - Singleton
0059
0060   fMyCellParameterisation = CellParameterisation::Instance();
0061
0062   // Material : 1 is cytoplasm, 2 is nucleus
0063
0064   G4int matVoxelPRE = -1;
0065   G4int matVoxelPOST = -1;
0066
0067   const G4StepPoint* preStep = aStep->GetPreStepPoint();
0068   const G4StepPoint* postStep = aStep->GetPostStepPoint();
0069   const G4Track* track = aStep->GetTrack();
0070
0071   const G4LogicalVolume* preVolume =
0072     preStep->GetPhysicalVolume()->GetLogicalVolume();
0073
0074   const G4LogicalVolume* postVolume = nullptr;
0075   if(postStep->GetPhysicalVolume())
0076   {
0077     postVolume = postStep->GetPhysicalVolume()->GetLogicalVolume();
0078   }
0079   const G4ParticleDefinition* particle =
0080     track->GetDynamicParticle()->GetDefinition();
0081
0082   G4int preReplicaNumber = preStep->GetTouchableHandle()->GetReplicaNumber();
0083   G4double edep = aStep->GetTotalEnergyDeposit();
0084
0085   if (preReplicaNumber>0)
0086   {
0087     matVoxelPRE =  fMyCellParameterisation->GetTissueType(preReplicaNumber);
0088   }
0089
0090   if(postVolume)
0091   {
0092     G4int postReplicaNumber = postStep->GetTouchableHandle()->GetReplicaNumber();
0093     if (postReplicaNumber>0)
0094     {
0095       matVoxelPOST = fMyCellParameterisation->GetTissueType(postReplicaNumber);
0096     }
0097   }
0098
0099   // COUNT GAS DETECTOR HITS
0100
0101   if (particle == G4Alpha::AlphaDefinition())
0102   {
0103     if(postVolume == fDetector->GetLogicalIsobutane() &&
0104        ((preVolume == fDetector->GetLogicalCollDetYoke())
0105     ||
0106     (preVolume == fDetector->GetLogicalCollDetGap4())
0107     ||
0108     (preVolume == fDetector->GetLogicalCollDetGap4())))
0109     {
0110       fRun->AddNbOfHitsGas();
0111     }
0112
0113     // STOPPING POWER AND BEAM SPOT SIZE AT CELL ENTRANCE
0114     if(preVolume == fDetector->GetLogicalPolyprop() &&
0115        ( (postVolume == fDetector->GetLogicalKgm()) ||
0116      (matVoxelPOST == 1)) )
0117     {
0118       G4double deltaE = preStep->GetKineticEnergy()
0119     - postStep->GetKineticEnergy();
0120       if(deltaE > 0.0)
0121       {
0122     //Fill ntupleid=1
0123     man->FillNtupleDColumn(1,0,preStep->GetKineticEnergy()/keV);
0124     man->FillNtupleDColumn(1,1,deltaE*micrometer/(keV*aStep->GetStepLength()));
0125     man->AddNtupleRow(1);
0126       }
0127
0128       // Average dE over step suggested by Michel Maire
0129       G4ThreeVector coord1 = preStep->GetPosition();
0130       const G4AffineTransform transformation1 =
0131     preStep->GetTouchable()->GetHistory()->GetTopTransform();
0132       G4ThreeVector localPosition1 = transformation1.TransformPoint(coord1);
0133
0134       G4ThreeVector coord2 = postStep->GetPosition();
0135       const G4AffineTransform transformation2 =
0136     postStep->GetTouchable()->GetHistory()->GetTopTransform();
0137       G4ThreeVector localPosition2 = transformation2.TransformPoint(coord2);
0138
0139       G4ThreeVector localPosition =
0140     localPosition1 + G4UniformRand()*(localPosition2-localPosition1);
0141
0142       //Fill ntupleid=2
0143       man->FillNtupleDColumn(2,0,localPosition.x()/micrometer);
0144       man->FillNtupleDColumn(2,1,localPosition.y()/micrometer);
0145       man->AddNtupleRow(2);
0146     }
0147
0148     // ALPHA RANGE
0149     if (postStep->GetKineticEnergy() < eV &&
0150     ( (matVoxelPOST==1) ||
0151       (postVolume == fDetector->GetLogicalKgm()) ||
0152       (matVoxelPOST==2) ) )
0153     {
0154       //Fill ntupleid=3
0155       man->FillNtupleDColumn(3,0,postStep->GetPosition().x()/micrometer);
0156       man->FillNtupleDColumn(3,1,postStep->GetPosition().y()/micrometer);
0157       man->FillNtupleDColumn(3,2,postStep->GetPosition().z()/micrometer);
0158       man->AddNtupleRow(3);
0159     }
0160
0161     // TOTAL DOSE DEPOSIT AND DOSE DEPOSIT WITHIN A PHANTOM VOXEL
0162     // FOR ALL PARTICLES
0163   }
0164
0165   if (matVoxelPRE  == 2)
0166   {
0167     G4double dose = (edep/joule)/(fRun->GetMassNucleus()/kg);
0168     fRun->AddDoseN(dose);
0169     fRun->AddDoseBox(preReplicaNumber, edep/eV);
0170   }
0171   else if (matVoxelPRE  == 1)
0172   {
0173     G4double dose = (edep/joule)/(fRun->GetMassCytoplasm()/kg);
0174     fRun->AddDoseC(dose);
0175     fRun->AddDoseBox(preReplicaNumber, edep/eV);
0176   }
0177
0178   // PROTECTION AGAINST POSSIBLE MSC LOOPS FOR e-
0179
0180   //  if ( edep/MeV<1e-25 && particle == G4Electron::Electron())
0181   //  {
0182       //aStep->GetTrack()->SetTrackStatus(fStopAndKill);
0183       /*
0184     G4cout << "*** Warning *** : msc loop for "
0185     << track->GetDefinition()->GetParticleName()
0186     << " in " <<
0187     postPoint->GetTouchableHandle()->GetVolume()->GetName() << G4endl;
0188       */
0189   //  }
0190 }