Warning, file /geant4/examples/advanced/dna/dsbandrepair/src/ChemGeoImport.cc was not indexed
or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030 #include "ChemGeoImport.hh"
0031 #include "G4Filesystem.hh"
0032 #include "G4DNAMolecule.hh"
0033
0034 ChemGeoImport::ChemGeoImport()
0035 {
0036 GetVoxelDefFilePathList();
0037 fpGun = new UserMoleculeGun();
0038 }
0039
0040
0041
0042 ChemGeoImport::~ChemGeoImport()
0043 {
0044 if(fpGun)
0045 delete fpGun;
0046 }
0047
0048
0049
0050 void ChemGeoImport::InsertMoleculeInWorld()
0051 {
0052
0053
0054 if(fIsParsed)
0055 {
0056
0057
0058
0059 for(G4int i=0, ie=fMolecules.size(); i<ie; i++)
0060 {
0061
0062
0063 G4String name = fMolecules[i].fName;
0064
0065 G4ThreeVector moleculePosition = fMolecules[i].fPosition;
0066 G4int copyNum = fMolecules[i].fCopyNumber;
0067
0068 G4int strand = fMolecules[i].fStrand;
0069 ChemMolecule Bmolecule(name, copyNum, moleculePosition, strand, -1, -1, -1);
0070
0071 if(name=="phosphate1" || name=="phosphate2")
0072 {
0073 name=G4Phosphate::Definition()->GetName();
0074 Bmolecule.fName=name;
0075 }
0076 else if(name=="deoxyribose1" || name=="deoxyribose2")
0077 {
0078 name=G4Deoxyribose::Definition()->GetName();
0079 Bmolecule.fName=name;
0080 }
0081 else if(name=="base_adenine")
0082 {
0083 name=G4Adenine::Definition()->GetName();
0084 Bmolecule.fName=name;
0085 }
0086 else if(name=="base_guanine")
0087 {
0088 name=G4Guanine::Definition()->GetName();
0089 Bmolecule.fName=name;
0090 }
0091 else if(name=="base_thymine")
0092 {
0093 name=G4Thymine::Definition()->GetName();
0094 Bmolecule.fName=name;
0095 }
0096 else if(name=="base_cytosine")
0097 {
0098 name=G4Cytosine::Definition()->GetName();
0099 Bmolecule.fName=name;
0100 }
0101 else if(name=="histone")
0102 {
0103 name=G4Histone::Definition()->GetName();
0104 Bmolecule.fName=name;
0105 }
0106 else if(name=="solvatedElectron")
0107 {
0108 name=G4Electron_aq::Definition()->GetName();
0109 }
0110 else if(name=="water")
0111 {
0112 name=G4H2O::Definition()->GetName();
0113 }
0114 else
0115 {
0116 G4String msg =
0117 "The name "+ name+" is not specified in the listed chemical molecules";
0118 G4Exception("ChemGeoImport::BuildGeometry", "", FatalException, msg);
0119 }
0120
0121
0122 G4bool toBeRemoved = IsMoleculeInTheRemoveTable(Bmolecule);
0123 if(!toBeRemoved)
0124 {
0125
0126
0127 if(name != G4H2O::Definition()->GetName() )
0128 fpGun->AddMolecule(name, moleculePosition, 1.e-12*s, copyNum, strand);
0129 else
0130 fpGun->AddWaterMolecule(moleculePosition, fMolecules.at(i).fTrackId,
0131 ElectronicModification(fMolecules.at(i).fState),
0132 fMolecules.at(i).fElectronicLevel);
0133 }
0134
0135 }
0136 G4DNAChemistryManager::Instance()->SetGun(fpGun);
0137 }
0138 else
0139 {
0140 G4String msg =
0141 "ChemGeoImport::InsertMoleculeInWorld: The parse method needs to be called first.";
0142 G4Exception("ChemGeoImport::ChemGeoImport::InsertMoleculeInWorld", "",FatalException, msg);
0143 }
0144 }
0145
0146
0147
0148 void ChemGeoImport::Reset()
0149 {
0150
0151 if(fpGun){
0152 delete fpGun;
0153 fpGun = new UserMoleculeGun();
0154 }
0155 fMolecules.clear();
0156 fMolecules.shrink_to_fit();
0157 fToBeRemovedMol.clear();
0158 fIsParsed = false;
0159 fFactor = 1.;
0160 }
0161
0162
0163
0164 void ChemGeoImport::ParseFiles(const G4String& chemInputFile)
0165 {
0166 G4fs::path aP{std::string(chemInputFile)};
0167 if (G4fs::exists(aP)) {
0168 Reset();
0169 ParseChemInputFile(chemInputFile);
0170 auto geoPathFileName = GetVoxelDefFilePath(fGeoNameFromChemInput);
0171
0172 ParseGeoFile(geoPathFileName);
0173 fIsParsed = true;
0174 }
0175 }
0176
0177
0178
0179 void ChemGeoImport::ParseChemInputFile(const G4String& fileName)
0180 {
0181
0182 std::ifstream file;
0183 file.open(fileName.c_str() );
0184
0185 if(!file.good() )
0186 {
0187
0188 G4String msg = fileName+" could not be opened";
0189 G4Exception("ChemGeoImport::ParseChemInputFile", "", FatalException, msg);
0190 }
0191
0192
0193 G4String line;
0194
0195
0196 while(std::getline(file, line) )
0197 {
0198
0199 if(line.empty() )
0200 continue;
0201
0202
0203 std::istringstream issLine(line);
0204
0205
0206 G4String firstItem;
0207
0208
0209 issLine >> firstItem;
0210
0211
0212 if(firstItem=="#")
0213 continue;
0214
0215 else if(firstItem=="_input")
0216 {
0217 G4int type(-1), state(-1), electronicLevel(-1), parentTrackId(-1);
0218 G4double x, y, z;
0219 issLine >> type >> state >> electronicLevel;
0220 issLine >> x >> y >> z;
0221 issLine >> parentTrackId;
0222
0223 x *= fFactor*nm;
0224 y *= fFactor*nm;
0225 z *= fFactor*nm;
0226
0227 G4String name;
0228 if(type==1)
0229 name="water";
0230 else if(type==2)
0231 name="solvatedElectron";
0232 else
0233 {
0234 G4ExceptionDescription description;
0235 description << "The type " << type <<" is not recognized";
0236 G4Exception("ChemGeoImport::ParseFile", "Fatal", FatalException, description, "");
0237 }
0238
0239 ChemMolecule molecule(name, -1, G4ThreeVector(x,y,z), -1,
0240 state, electronicLevel, parentTrackId);
0241
0242 fMolecules.push_back(molecule);
0243 }
0244
0245 else if(firstItem=="_remove")
0246 {
0247 G4String name;
0248 issLine >> name;
0249
0250 G4int copyNumber;
0251 issLine >> copyNumber;
0252
0253 G4int strand;
0254 issLine >> strand;
0255
0256 fToBeRemovedMol.push_back(ChemMolecule(name,copyNumber,G4ThreeVector(),strand,-1,-1,-1));
0257 }
0258
0259 else if(firstItem=="_eventNum")
0260 {
0261
0262 }
0263
0264 else if(firstItem=="_voxelType")
0265 {
0266 issLine >> fGeoNameFromChemInput;
0267 }
0268
0269 else if(firstItem=="_voxelCopyNumber")
0270 {
0271
0272 }
0273
0274 else if(firstItem=="_Version")
0275 {
0276
0277 }
0278
0279 else
0280 {
0281
0282 G4String msg =
0283 firstItem+" is not defined in the parser. Check the input file: "+fileName+".";
0284 G4Exception("ChemGeoImport::ParseChemInputFile", "Geo_WrongParse",FatalException, msg);
0285 }
0286 }
0287 file.close();
0288 }
0289
0290
0291
0292 void ChemGeoImport::ParseGeoFile(const G4String& fileName)
0293 {
0294
0295 std::ifstream file(fileName.c_str());
0296
0297
0298 if(!file.is_open() )
0299 {
0300
0301 G4String msg = fileName+" could not be opened";
0302 G4Exception("ChemGeoImport::ParseGeoFile", "", FatalException, msg);
0303 }
0304
0305
0306 G4String line;
0307
0308
0309 while(std::getline(file, line) )
0310 {
0311
0312 if(line.empty() )
0313 continue;
0314
0315
0316 std::istringstream issLine(line);
0317
0318
0319 G4String firstItem;
0320
0321
0322 issLine >> firstItem;
0323
0324
0325 if(firstItem=="#")
0326 continue;
0327
0328
0329 else if(firstItem=="_Name")
0330 {
0331 G4String name;
0332 issLine >> name;
0333 }
0334 else if(firstItem=="_Size")
0335 {
0336 G4double size;
0337 issLine >> size;
0338 size *= fFactor*nm;
0339
0340 fSize = size;
0341 }
0342 else if(firstItem=="_Number")
0343 {
0344
0345 }
0346 else if(firstItem=="_Radius")
0347 {
0348
0349 }
0350 else if(firstItem=="_Version")
0351 {
0352
0353 }
0354 else if(firstItem=="_pl")
0355 {
0356 G4String name;
0357 issLine >> name;
0358
0359 G4String material;
0360 issLine >> material;
0361
0362 G4int strand;
0363 issLine >> strand;
0364
0365 G4int copyNumber;
0366 issLine >> copyNumber;
0367
0368 G4double x;
0369 issLine >> x;
0370 x *= fFactor*nm;
0371
0372 G4double y;
0373 issLine >> y;
0374 y *= fFactor*nm;
0375
0376 G4double z;
0377 issLine >> z;
0378 z *= fFactor*nm;
0379
0380 ChemMolecule molecule(name, copyNumber, G4ThreeVector(x, y, z), strand, -1, -1, -1);
0381
0382 fMolecules.push_back(molecule);
0383 }
0384
0385 else
0386 {
0387
0388 G4String msg =
0389 firstItem+" is not defined in the parser. Check the input file: "+fileName+".";
0390 G4Exception("ChemGeoImport::ParseGeoFile", "Geo_WrongParse", FatalException, msg);
0391 }
0392 }
0393 file.close();
0394 }
0395
0396
0397
0398 G4bool ChemGeoImport::IsMoleculeInTheRemoveTable(const ChemMolecule& molecule)
0399 {
0400 if(std::find(fToBeRemovedMol.begin(),fToBeRemovedMol.end(),molecule) != fToBeRemovedMol.end())
0401 return true;
0402 else
0403 return false;
0404 }
0405
0406
0407
0408 G4String ChemGeoImport::GetVoxelDefFilePath(G4String bareName)
0409 {
0410 G4String strRes = "";
0411 for (auto const &entry : fVoxelDefFilesList) {
0412 G4fs::path voxelP{std::string(entry)};
0413 if (voxelP.stem().string() == bareName) {
0414 strRes = entry;
0415 }
0416 }
0417 return strRes;
0418 }
0419
0420
0421
0422 void ChemGeoImport::GetVoxelDefFilePathList()
0423 {
0424 G4fs::path thisP = G4fs::current_path();
0425 G4bool doesWantedFileExist = false;
0426 for (const auto &entry : G4fs::directory_iterator(thisP)){
0427 if (entry.path().filename() == "imp.info") {
0428 std::ifstream file(entry.path().c_str());
0429 if(!file.good() ){
0430 G4String msg =
0431 "File imp.info is broken. Check its content or try to rerun the PhysicalStage?";
0432 G4Exception("ChemGeoImport::GetVoxelDefFilePathList()", "", FatalException, msg);
0433 }
0434 doesWantedFileExist = true;
0435 G4String line;
0436 while(std::getline(file, line) ){
0437 std::istringstream iss(line);
0438 G4String flag;
0439 G4String voxelDefFile;
0440 iss >> flag;
0441 if ( flag == "_geovolxelpath") {
0442 iss >> voxelDefFile;
0443 fVoxelDefFilesList.insert(voxelDefFile);
0444 }
0445 }
0446 file.close();
0447 }
0448 }
0449
0450 if (!doesWantedFileExist) {
0451 G4String msg = "File imp.info does not exist. Did you run the Physical Stage?";
0452 G4Exception("ChemGeoImport::GetVoxelDefFilePathList()", "", FatalException, msg);
0453 }
0454 }
0455
0456