Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-02-23 09:22:01

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 PrimaryGeneratorSourceGRASCSV.cc
0028 /// \brief Primary Generator source class implementation for GRAS CSV phase space for Molecular DNA
0029 /// simulation
0030 
0031 #include "PrimaryGeneratorSourceGRASCSV.hh"
0032 
0033 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0034 
0035 PrimaryGeneratorSourceGRASCSV::PrimaryGeneratorSourceGRASCSV(const G4String& filename)
0036 {
0037   // set a default buffer size
0038   fEndOfFile = false;
0039 
0040   // open file
0041   fInputFile.open(filename);
0042   if (fInputFile.is_open()) {
0043     G4cout << "*** Opening particle source file " << filename << " ***" << G4endl;
0044     // find the starting particle source box
0045     G4String line = "";
0046     G4int totalevents = 0;
0047     while (std::getline(fInputFile, line)) {
0048       if (line.find("TOTAL_EVENTS") != std::string::npos) {
0049         // Create a stringstream of the current line
0050         std::stringstream ss(line);
0051         G4String tmp = "";
0052         ss >> tmp >> tmp >> totalevents;
0053       }
0054       if (line.find("TWO-STAGE PARTICLE PHASE SPACE OUTPUT FILE: PARTICLE PHASE-SPACE INFORMATION")
0055           != std::string::npos)
0056       {
0057         for (G4int i = 0; i < 18; i++) {
0058           std::getline(fInputFile, line);
0059           // Print the header information
0060           G4cout << line << G4endl;
0061         }
0062         break;
0063       }
0064     }
0065     // set the number of particles
0066     fnParticles = totalevents;
0067   }
0068   else {
0069     G4cout << "*** Warning: can't open particle source file in reading mode ***" << G4endl;
0070   }
0071 }
0072 
0073 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0074 
0075 PrimaryGeneratorSourceGRASCSV::~PrimaryGeneratorSourceGRASCSV()
0076 {
0077   // close file and clear structures
0078   if (fInputFile.is_open()) {
0079     fInputFile.close();
0080   }
0081   fPrimaryList.clear();
0082 }
0083 
0084 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0085 
0086 G4double PrimaryGeneratorSourceGRASCSV::RecomputeNParticles()
0087 {
0088   // to be implemented
0089   // the variable in CSV has to be recomputed reading all the file
0090   return fnParticles;
0091 }
0092 
0093 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0094 
0095 Primary* PrimaryGeneratorSourceGRASCSV::GetPrimary()
0096 {
0097   if (fInputFile.peek() == EOF || fEndOfFile == true) return nullptr;
0098   // Read in primaries if the primary list is empty
0099   if (fPrimaryList.size() == 0) {
0100     // Read 100 primaries at a time
0101     G4int num = fBufferSize;
0102     if (num > fnParticles) {
0103       num = fnParticles;
0104     }
0105 
0106     for (G4int i = 0; i < num; i++) {
0107       if (fInputFile.peek() == EOF) {
0108         // End of file
0109         fEndOfFile = true;
0110         break;
0111       }
0112       G4String line = "";
0113       std::getline(fInputFile, line);
0114       if (line.find("TWO-STAGE PARTICLE PHASE SPACE OUTPUT FILE: PARTICLE PHASE-SPACE INFORMATION")
0115           != std::string::npos)
0116       {
0117         const G4int lines = 18;
0118         for (G4int j = 0; j < lines; j++) {
0119           std::getline(fInputFile, line);
0120         }
0121         i--;
0122         continue;
0123       }
0124       if (line.find("Block") != std::string::npos || line.find("File") != std::string::npos
0125           || line.find("End") != std::string::npos || line.find("*") != std::string::npos)
0126       {
0127         // End of Block, skip
0128         i--;
0129         continue;
0130       }
0131       std::stringstream ss(line);
0132       std::vector<G4String> tempVect;
0133       G4String tempString = "";
0134       while (ss >> tempString)
0135         tempVect.push_back(tempString);
0136 
0137       try {
0138         // Define position and momentum direction vectors
0139         G4ThreeVector pos = G4ThreeVector(std::stod(tempVect.at(10)), std::stod(tempVect.at(11)),
0140                                           std::stod(tempVect.at(12)));
0141         G4ThreeVector momDir = G4ThreeVector(std::stod(tempVect.at(16)), std::stod(tempVect.at(17)),
0142                                              std::stod(tempVect.at(18)));
0143         G4int PDGEncoding = std::stoi(tempVect.at(2));
0144         G4double KE = 0;
0145         KE = std::stod(tempVect.at(8));
0146 
0147         // Generate primary particle
0148         auto* primary = new Primary();
0149         primary->SetName(PDGEncoding);
0150         primary->SetPosition(pos);
0151         primary->SetMomentumDirection(momDir);
0152         primary->SetEnergy(KE);
0153 
0154         // Populate the primaries list
0155         fPrimaryList.push_back(primary);
0156       }
0157       catch (const std::invalid_argument&) {
0158         G4cerr << "Phase Space reading error: invalid_argument" << G4endl;
0159       }
0160       catch (const std::out_of_range&) {
0161         G4cerr << "Phase Space reading error: out_of_range" << G4endl;
0162       }
0163     }
0164   }
0165 
0166   // Get first element and delete it so it's not reused
0167   Primary* primary = nullptr;
0168   if (fPrimaryList.size() > 0) {
0169     primary = fPrimaryList.front();
0170     fPrimaryList.pop_front();
0171   }
0172   return primary;
0173 }