dsbandrepair/src/PhysChemIO.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 //
0027 /// \file PhysChemIO.cc
0028 /// \brief Implementation of the PhysChemIO class
0029
0030 #include "PhysChemIO.hh"
0031 #include "SteppingAction.hh"
0032 #include "Analysis.hh"
0033 #include "G4Track.hh"
0034 #include "G4NavigationHistory.hh"
0035 #include "G4RunManager.hh"
0036
0037 #ifdef USE_MPI
0038 #include "G4MPImanager.hh"
0039 #endif
0040
0041 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0042
0043 PhysChemIO::PhysChemIO(SteppingAction* stepAction) : G4VPhysChemIO(),
0044 fSteppingAction(stepAction)
0045 {;}
0046
0047 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0048
0049 void PhysChemIO::CreateWaterMolecule(G4int electronicModif, G4int electronicLevel,
0050                                     G4double /*energy*/,
0051                                     const G4Track* theIncomingTrack)
0052 {
0053     //L.T. Anh:  to correct electronicLevel in G4DNAChemistryManager,
0054     //see in G4DNAChemistryManager::CreateWaterMolecule
0055     electronicLevel = 4 - electronicLevel;
0056     //Rel pos
0057     G4ThreeVector relPos;
0058     auto touchable = theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable();
0059     relPos = touchable->GetHistory()->GetTopTransform().TransformPoint(theIncomingTrack->GetPosition());
0060
0061     // Get the flag of the current volume
0062     G4int volumeFlag =fSteppingAction->SetupVolumeFlag(
0063         theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetVolume()->GetName());
0064
0065     if(    volumeFlag == 161  // voxelStraight
0066         || volumeFlag == 162 // voxelRight
0067         || volumeFlag == 163 // voxelLeft
0068         || volumeFlag == 164 // voxelUp
0069         || volumeFlag == 165 // voxelDown
0070         || volumeFlag == 261  // voxelStraight2
0071         || volumeFlag == 262 // voxelRight2
0072         || volumeFlag == 263 // voxelLeft2
0073         || volumeFlag == 264 // voxelUp2
0074         || volumeFlag == 265) // voxelDown2
0075     {
0076         // Get the volume copy number
0077         G4int volumeCpNum = touchable->GetCopyNumber();
0078         //theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetVolume()->GetUserID();
0079
0080         // Default flag values
0081         G4String motherVolumeName = "";
0082         G4int motherVolumeFlag = -1;
0083         G4int motherVolumeCpNum = -1;
0084
0085         // Mother volume informations
0086
0087         // Be sure there is a mother volume to ask for
0088         if(theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetHistoryDepth() >0)
0089         {
0090             G4VPhysicalVolume* motherVol = theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetVolume(1);
0091
0092             // General infos
0093             motherVolumeName = motherVol->GetName();
0094             motherVolumeFlag = fSteppingAction->SetupVolumeFlag(motherVolumeName);
0095             motherVolumeCpNum = motherVol->GetCopyNo();
0096         }
0097         G4int eventId = G4RunManager::GetRunManager()->GetCurrentEvent()->GetEventID();
0098 #ifdef USE_MPI
0099         auto g4MPI = G4MPImanager::GetManager();
0100         if (g4MPI->IsSlave()) { // update eventID only for slave, cause rank_master=0
0101             G4int rank = g4MPI->GetRank();
0102             eventId += g4MPI->GetEventsInMaster() + (rank-1)*g4MPI->GetEventsInSlave();
0103         }
0104 #endif
0105
0106         InfoForChemGeo aInfo;
0107         aInfo.fType = 1; // water=1
0108         aInfo.fState = electronicModif ;
0109         aInfo.fElectronicLevel = electronicLevel ;
0110         aInfo.fX = theIncomingTrack->GetPosition().x()/nm;
0111         aInfo.fY = theIncomingTrack->GetPosition().y()/nm;
0112         aInfo.fZ = theIncomingTrack->GetPosition().z()/nm;
0113         aInfo.fParentTrackID = theIncomingTrack->GetTrackID() ;
0114         aInfo.fEventNumber = eventId;
0115         aInfo.fVolume = volumeFlag ;
0116         aInfo.fVolumeCopyNumber = volumeCpNum;
0117         aInfo.fMotherVolume = motherVolumeFlag ;
0118         aInfo.fMotherVolumeCopyNumber = motherVolumeCpNum ;
0119         aInfo.fRelX = relPos.x()/nm;
0120         aInfo.fRelY = relPos.y()/nm;
0121         aInfo.fRelZ = relPos.z()/nm;
0122         Analysis::GetAnalysis()->AddInfoForChemGeo(aInfo);
0123     }
0124 }
0125
0126 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0127
0128 void PhysChemIO::CreateSolvatedElectron(const G4Track* theIncomingTrack, G4ThreeVector* finalPosition)
0129 {
0130     G4ThreeVector pos;
0131     if(finalPosition) pos = *finalPosition;
0132     else pos = theIncomingTrack->GetPosition();
0133
0134     // Rel pos
0135     G4ThreeVector relPos;
0136     const G4VTouchable* touchable = theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable();
0137     relPos = touchable->GetHistory()->GetTopTransform().TransformPoint(pos);
0138
0139     // Current volume infos
0140     G4int volumeFlag = fSteppingAction->SetupVolumeFlag(
0141         theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetVolume()->GetName());
0142     if(    volumeFlag == 161  // voxelStraight
0143         || volumeFlag == 162 // voxelRight
0144         || volumeFlag == 163 // voxelLeft
0145         || volumeFlag == 164 // voxelUp
0146         || volumeFlag == 165 // voxelDown
0147         || volumeFlag == 261  // voxelStraight2
0148         || volumeFlag == 262 // voxelRight2
0149         || volumeFlag == 263 // voxelLeft2
0150         || volumeFlag == 264 // voxelUp2
0151         || volumeFlag == 265) // voxelDown2
0152     {
0153         G4int volumeCpNum = touchable->GetCopyNumber();
0154
0155         G4String motherVolumeName = "";
0156         G4int motherVolumeFlag = -1;
0157         G4int motherVolumeCpNum = -1;
0158
0159         // Mother volume informations
0160         // Be sure there is a mother volume to ask for
0161         if(theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetHistoryDepth() >0)
0162         {
0163             G4VPhysicalVolume* motherVol = theIncomingTrack->GetStep()->GetPreStepPoint()->GetTouchable()->GetVolume(1);
0164             // General infos
0165             motherVolumeName = motherVol->GetName();
0166             motherVolumeFlag = fSteppingAction->SetupVolumeFlag(motherVolumeName);
0167             motherVolumeCpNum = motherVol->GetCopyNo();
0168         }
0169         G4int eventId = G4RunManager::GetRunManager()->GetCurrentEvent()->GetEventID();
0170 #ifdef USE_MPI
0171         auto g4MPI = G4MPImanager::GetManager();
0172         if (g4MPI->IsSlave()) { // update eventID only for slave, cause rank_master=0
0173             G4int rank = g4MPI->GetRank();
0174             eventId += g4MPI->GetEventsInMaster() + (rank-1)*g4MPI->GetEventsInSlave();
0175         }
0176 #endif
0177
0178         InfoForChemGeo aInfo;
0179         aInfo.fType = 2; // / solvated electron=2
0180         aInfo.fState = -1; // no state for solvated electron
0181         aInfo.fElectronicLevel = -1; // no electronic level for solvated electron
0182         aInfo.fX = pos.x()/nm;
0183         aInfo.fY = pos.y()/nm;
0184         aInfo.fZ = pos.z()/nm;
0185         aInfo.fParentTrackID = theIncomingTrack->GetTrackID() ;
0186         aInfo.fEventNumber = eventId;
0187         aInfo.fVolume = volumeFlag ;
0188         aInfo.fVolumeCopyNumber = volumeCpNum ;
0189         aInfo.fMotherVolume = motherVolumeFlag ;
0190         aInfo.fMotherVolumeCopyNumber =  motherVolumeCpNum ;
0191         aInfo.fRelX = relPos.x()/nm;
0192         aInfo.fRelY = relPos.y()/nm;
0193         aInfo.fRelZ = relPos.z()/nm;
0194         Analysis::GetAnalysis()->AddInfoForChemGeo(aInfo);
0195     }
0196 }
0197
0198 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......