Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 /* 26 /* 27 * G4ReactionTableMessenger.cc 27 * G4ReactionTableMessenger.cc 28 * 28 * 29 * Created on: Sep 14, 2015 29 * Created on: Sep 14, 2015 30 * Author: mkaramit 30 * Author: mkaramit 31 */ 31 */ 32 32 33 #include <G4ReactionTableMessenger.hh> 33 #include <G4ReactionTableMessenger.hh> 34 #include <G4DNAMolecularReactionTable.hh> 34 #include <G4DNAMolecularReactionTable.hh> 35 #include <G4UIcmdWithAString.hh> 35 #include <G4UIcmdWithAString.hh> 36 #include <G4UIcmdWithADoubleAndUnit.hh> 36 #include <G4UIcmdWithADoubleAndUnit.hh> 37 #include <G4UnitsTable.hh> 37 #include <G4UnitsTable.hh> 38 #include <G4SystemOfUnits.hh> 38 #include <G4SystemOfUnits.hh> 39 #include <G4UIcmdWithoutParameter.hh> 39 #include <G4UIcmdWithoutParameter.hh> 40 #include <G4UIcmdWithABool.hh> 40 #include <G4UIcmdWithABool.hh> 41 41 42 //-------------------------------------------- 42 //------------------------------------------------------------------------------ 43 43 44 G4ReactionTableMessenger::G4ReactionTableMesse 44 G4ReactionTableMessenger::G4ReactionTableMessenger(G4DNAMolecularReactionTable* table) 45 : 45 : 46 fpTable(table) 46 fpTable(table) 47 , fpActivateReactionUI(new G4UIcmdWithoutPar 47 , fpActivateReactionUI(new G4UIcmdWithoutParameter("/chem/reaction/UI", this)) 48 { 48 { 49 fpNewDiffContReaction = new G4UIcmdWithAStri 49 fpNewDiffContReaction = new G4UIcmdWithAString("/chem/reaction/new", this); 50 fpAddReaction = new G4UIcmdWithAString("/che 50 fpAddReaction = new G4UIcmdWithAString("/chem/reaction/add", this); 51 fpPrintTable = new G4UIcmdWithoutParameter(" 51 fpPrintTable = new G4UIcmdWithoutParameter("/chem/reaction/print", this); 52 } 52 } 53 53 54 //-------------------------------------------- 54 //------------------------------------------------------------------------------ 55 55 56 G4ReactionTableMessenger::~G4ReactionTableMess 56 G4ReactionTableMessenger::~G4ReactionTableMessenger() 57 { 57 { 58 delete fpNewDiffContReaction; 58 delete fpNewDiffContReaction; 59 delete fpAddReaction; 59 delete fpAddReaction; 60 delete fpPrintTable; 60 delete fpPrintTable; 61 } 61 } 62 62 63 //-------------------------------------------- 63 //------------------------------------------------------------------------------ 64 void G4ReactionTableMessenger::SetNewValue(G4U 64 void G4ReactionTableMessenger::SetNewValue(G4UIcommand* command, 65 G4S 65 G4String newValue) 66 { 66 { 67 if(command == fpActivateReactionUI.get()) 67 if(command == fpActivateReactionUI.get()) 68 { 68 { 69 //assert(false); 69 //assert(false); 70 fpTable->Reset();//release reaction data 70 fpTable->Reset();//release reaction data 71 } 71 } 72 72 73 if(command == fpNewDiffContReaction) 73 if(command == fpNewDiffContReaction) 74 { 74 { 75 std::istringstream iss(newValue); 75 std::istringstream iss(newValue); 76 76 77 G4String species1; 77 G4String species1; 78 iss >> species1; 78 iss >> species1; 79 79 80 G4String species2; 80 G4String species2; 81 iss >> species2; 81 iss >> species2; 82 82 83 double reactionRate; 83 double reactionRate; 84 iss >> reactionRate; 84 iss >> reactionRate; 85 85 86 // G4String reactionRateUnit; 86 // G4String reactionRateUnit; 87 // iss >> reactionRateUnit; 87 // iss >> reactionRateUnit; 88 88 89 double dimensionedReactionRate = reactionR 89 double dimensionedReactionRate = reactionRate * (1e-3 * m3 / (mole * s)); 90 // G4UIcmdWithADoubleAndUnit::ConvertTo 90 // G4UIcmdWithADoubleAndUnit::ConvertToDimensionedDouble((reactionRate 91 // + G4String(" ") + reactionRateUn 91 // + G4String(" ") + reactionRateUnit).c_str()); 92 92 93 auto reactionData = 93 auto reactionData = 94 new G4DNAMolecularReactionData(dimensi 94 new G4DNAMolecularReactionData(dimensionedReactionRate, 95 species 95 species1, 96 species 96 species2); 97 97 98 // G4String productionRate; 98 // G4String productionRate; 99 // iss >> productionRate; 99 // iss >> productionRate; 100 // 100 // 101 // if(productionRate != "" && productionRat 101 // if(productionRate != "" && productionRate != "X") 102 // { 102 // { 103 // double prodRateReal = G4UIcommand::Con 103 // double prodRateReal = G4UIcommand::ConvertToDouble(productionRate); 104 // 104 // 105 // if(prodRateReal == 0 || isnan(prodRate 105 // if(prodRateReal == 0 || isnan(prodRateReal)) 106 // { 106 // { 107 // G4Exception("G4ReactionTableMessenge 107 // G4Exception("G4ReactionTableMessenger", 108 // "WRONG_PRODUCTION_RATE", 108 // "WRONG_PRODUCTION_RATE", 109 // FatalException, ""); 109 // FatalException, ""); 110 // } 110 // } 111 // 111 // 112 // double dimensionedProductionRate = pro 112 // double dimensionedProductionRate = prodRateReal 113 // * (1e-3 * m3 / (mole * s)); 113 // * (1e-3 * m3 / (mole * s)); 114 // reactionData->SetProductionRate(dimens 114 // reactionData->SetProductionRate(dimensionedProductionRate); 115 // } 115 // } 116 116 117 while(!iss.eof()) 117 while(!iss.eof()) 118 { 118 { 119 G4String product; 119 G4String product; 120 iss >> product; 120 iss >> product; 121 121 122 if(!product.empty()) 122 if(!product.empty()) 123 { 123 { 124 reactionData->AddProduct(product); 124 reactionData->AddProduct(product); 125 } 125 } 126 else 126 else 127 { 127 { 128 break; 128 break; 129 } 129 } 130 }; 130 }; 131 131 132 fpTable->SetReaction(reactionData); 132 fpTable->SetReaction(reactionData); 133 } 133 } 134 // else if(command == fpNewPartDiffContReacti 134 // else if(command == fpNewPartDiffContReactionByRadius) 135 // { 135 // { 136 // std::istringstream iss(newValue); 136 // std::istringstream iss(newValue); 137 // 137 // 138 // G4String species1; 138 // G4String species1; 139 // iss >> species1; 139 // iss >> species1; 140 // 140 // 141 // G4String species2; 141 // G4String species2; 142 // iss >> species2; 142 // iss >> species2; 143 // 143 // 144 // double reactionRate; 144 // double reactionRate; 145 // iss >> reactionRate; 145 // iss >> reactionRate; 146 // 146 // 147 //// G4String reactionRateUnit; 147 //// G4String reactionRateUnit; 148 //// iss >> reactionRateUnit; 148 //// iss >> reactionRateUnit; 149 // 149 // 150 //// G4String reactionRateUnit; 150 //// G4String reactionRateUnit; 151 //// iss >> reactionRateUnit; 151 //// iss >> reactionRateUnit; 152 // 152 // 153 // double reactionRadius; 153 // double reactionRadius; 154 // iss >> reactionRadius; 154 // iss >> reactionRadius; 155 // 155 // 156 //// G4String reactionRadiusUnit; 156 //// G4String reactionRadiusUnit; 157 //// iss >> reactionRadiusUnit; 157 //// iss >> reactionRadiusUnit; 158 // 158 // 159 // double dimensionedReactionRate = reactio 159 // double dimensionedReactionRate = reactionRate * (1e-3 * m3 / (mole * s)); 160 // 160 // 161 //// double dimensionedReactionRate = 161 //// double dimensionedReactionRate = 162 //// G4UIcmdWithADoubleAndUnit::Convert 162 //// G4UIcmdWithADoubleAndUnit::ConvertToDimensionedDouble((reactionRate 163 //// + " " + reactionRateUnit).c_st 163 //// + " " + reactionRateUnit).c_str()); 164 // 164 // 165 // double dimensionedReactionRadius = react 165 // double dimensionedReactionRadius = reactionRadius * nm; 166 //// G4UIcmdWithADoubleAndUnit::Convert 166 //// G4UIcmdWithADoubleAndUnit::ConvertToDimensionedDouble((reactionRadius 167 //// + " " + reactionRadiusUnit).c_ 167 //// + " " + reactionRadiusUnit).c_str()); 168 // 168 // 169 // G4DNAMolecularReactionData* reactionData 169 // G4DNAMolecularReactionData* reactionData = 170 // new G4DNAMolecularReactionData(dimen 170 // new G4DNAMolecularReactionData(dimensionedReactionRate, 171 // speci 171 // species1, 172 // speci 172 // species2); 173 // reactionData->SetPartiallyDiffusionContr 173 // reactionData->SetPartiallyDiffusionControlledReaction(dimensionedReactionRate, 174 // 174 // dimensionedReactionRadius); 175 // 175 // 176 // while(iss.eof() == false) 176 // while(iss.eof() == false) 177 // { 177 // { 178 // G4String product; 178 // G4String product; 179 // iss >> product; 179 // iss >> product; 180 // 180 // 181 // if(product != "") 181 // if(product != "") 182 // { 182 // { 183 // reactionData->AddProduct(product); 183 // reactionData->AddProduct(product); 184 // } 184 // } 185 // else 185 // else 186 // { 186 // { 187 // break; 187 // break; 188 // } 188 // } 189 // } 189 // } 190 // 190 // 191 // fpTable->SetReaction(reactionData); 191 // fpTable->SetReaction(reactionData); 192 // } 192 // } 193 // else if(command == fpNewPartDiffContReacti 193 // else if(command == fpNewPartDiffContReactionByReactionRate) 194 // { 194 // { 195 // std::istringstream iss(newValue); 195 // std::istringstream iss(newValue); 196 // 196 // 197 // G4String species1; 197 // G4String species1; 198 // iss >> species1; 198 // iss >> species1; 199 // 199 // 200 // G4String species2; 200 // G4String species2; 201 // iss >> species2; 201 // iss >> species2; 202 // 202 // 203 // double reactionRate; 203 // double reactionRate; 204 // iss >> reactionRate; 204 // iss >> reactionRate; 205 // 205 // 206 // // G4String reactionRateUnit; 206 // // G4String reactionRateUnit; 207 // // iss >> reactionRateUnit; 207 // // iss >> reactionRateUnit; 208 // 208 // 209 // // G4String reactionRateUnit; 209 // // G4String reactionRateUnit; 210 // // iss >> reactionRateUnit; 210 // // iss >> reactionRateUnit; 211 // 211 // 212 // double activationRate; 212 // double activationRate; 213 // iss >> activationRate; 213 // iss >> activationRate; 214 // 214 // 215 // // G4String reactionRadiusUnit; 215 // // G4String reactionRadiusUnit; 216 // // iss >> reactionRadiusUnit; 216 // // iss >> reactionRadiusUnit; 217 // 217 // 218 // double dimensionedReactionRate = reactio 218 // double dimensionedReactionRate = reactionRate * (1e-3 * m3 / (mole * s)); 219 // 219 // 220 // double dimensionedActivationRate = activ 220 // double dimensionedActivationRate = activationRate 221 // * (1e-3 * m3 / (mole * s)); 221 // * (1e-3 * m3 / (mole * s)); 222 // 222 // 223 // G4DNAMolecularReactionData* reactionData 223 // G4DNAMolecularReactionData* reactionData = 224 // new G4DNAMolecularReactionData(dimen 224 // new G4DNAMolecularReactionData(dimensionedReactionRate, 225 // speci 225 // species1, 226 // speci 226 // species2); 227 // reactionData-> 227 // reactionData-> 228 // SetPartiallyDiffusionControlledReactio 228 // SetPartiallyDiffusionControlledReactionByActivation(dimensionedReactionRate, 229 // 229 // dimensionedActivationRate); 230 // 230 // 231 // while(iss.eof() == false) 231 // while(iss.eof() == false) 232 // { 232 // { 233 // G4String product; 233 // G4String product; 234 // iss >> product; 234 // iss >> product; 235 // 235 // 236 // if(product != "") 236 // if(product != "") 237 // { 237 // { 238 // reactionData->AddProduct(product); 238 // reactionData->AddProduct(product); 239 // } 239 // } 240 // else 240 // else 241 // { 241 // { 242 // break; 242 // break; 243 // } 243 // } 244 // } 244 // } 245 // 245 // 246 // fpTable->SetReaction(reactionData); 246 // fpTable->SetReaction(reactionData); 247 // } 247 // } 248 else if(command == fpPrintTable) 248 else if(command == fpPrintTable) 249 { 249 { 250 fpTable->PrintTable(); 250 fpTable->PrintTable(); 251 } 251 } 252 else if(command == fpAddReaction) 252 else if(command == fpAddReaction) 253 { 253 { 254 std::istringstream iss(newValue); 254 std::istringstream iss(newValue); 255 255 256 //---------------------------------------- 256 //-------------------------------------------------------------------------- 257 // Reactants definition 257 // Reactants definition 258 258 259 G4String species1; 259 G4String species1; 260 iss >> species1; 260 iss >> species1; 261 261 262 G4String marker; 262 G4String marker; 263 iss >> marker; // peut etre +, ->, | 263 iss >> marker; // peut etre +, ->, | 264 264 265 G4String species2; 265 G4String species2; 266 266 267 if(marker == "+") 267 if(marker == "+") 268 { 268 { 269 iss >> species2; 269 iss >> species2; 270 iss >> marker; // peut etre ->, | 270 iss >> marker; // peut etre ->, | 271 } 271 } 272 272 273 //---------------------------------------- 273 //-------------------------------------------------------------------------- 274 274 275 auto reactionData = 275 auto reactionData = 276 new G4DNAMolecularReactionData 276 new G4DNAMolecularReactionData(0, 277 277 species1, 278 278 species2); 279 //fpTable->SetReaction(reactionData); 279 //fpTable->SetReaction(reactionData); 280 280 281 //---------------------------------------- 281 //-------------------------------------------------------------------------- 282 // Add products 282 // Add products 283 if(marker == "->") 283 if(marker == "->") 284 { 284 { 285 iss >> marker; // doit etre = species na 285 iss >> marker; // doit etre = species name 286 286 287 while(marker!="|" 287 while(marker!="|" 288 //&& marker!="" 288 //&& marker!="" 289 && !iss.eof() 289 && !iss.eof() 290 ) 290 ) 291 { 291 { 292 //G4cout << marker << G4endl; 292 //G4cout << marker << G4endl; 293 if(marker == "+") 293 if(marker == "+") 294 { 294 { 295 iss >> marker; // doit etre species 295 iss >> marker; // doit etre species name 296 continue; 296 continue; 297 } 297 } 298 if(marker != "H2O") 298 if(marker != "H2O") 299 { 299 { 300 reactionData->AddProduct(marker); 300 reactionData->AddProduct(marker); 301 } 301 } 302 302 303 iss >> marker; // peut etre species na 303 iss >> marker; // peut etre species name, +, | 304 }; 304 }; 305 } 305 } 306 306 307 // G4cout << "out of 1st loop" << G4endl; 307 // G4cout << "out of 1st loop" << G4endl; 308 308 309 //---------------------------------------- 309 //-------------------------------------------------------------------------- 310 // Add reaction rate method 310 // Add reaction rate method 311 G4String rateconst_method; 311 G4String rateconst_method; 312 iss >> rateconst_method; 312 iss >> rateconst_method; 313 if(rateconst_method == "Fix") 313 if(rateconst_method == "Fix") 314 { 314 { 315 iss >> marker; // must be | 315 iss >> marker; // must be | 316 double reactionRate; 316 double reactionRate; 317 iss >> reactionRate; 317 iss >> reactionRate; 318 318 319 double dimensionedReactionRate = reactio 319 double dimensionedReactionRate = reactionRate * (1e-3 * m3 / (mole * s)); 320 reactionData->SetObservedReactionRateCon 320 reactionData->SetObservedReactionRateConstant(dimensionedReactionRate); 321 reactionData->ComputeEffectiveRadius(); 321 reactionData->ComputeEffectiveRadius(); 322 G4String markerType; 322 G4String markerType; 323 iss >> markerType; // must be | 323 iss >> markerType; // must be | 324 if(markerType == "|") 324 if(markerType == "|") 325 { 325 { 326 G4int reactionType; 326 G4int reactionType; 327 iss >> reactionType; 327 iss >> reactionType; 328 if(reactionType == 1) 328 if(reactionType == 1) 329 { 329 { 330 reactionData->SetReactionType(reacti 330 reactionData->SetReactionType(reactionType); 331 } 331 } 332 } 332 } 333 333 334 334 335 // G4String productionRate; 335 // G4String productionRate; 336 // iss >> productionRate; 336 // iss >> productionRate; 337 // 337 // 338 // if(productionRate != "" && productionR 338 // if(productionRate != "" && productionRate != "X") 339 // { 339 // { 340 // double prodRateReal = G4UIcommand::C 340 // double prodRateReal = G4UIcommand::ConvertToDouble(productionRate); 341 // 341 // 342 // if(prodRateReal == 0 || isnan(prodRa 342 // if(prodRateReal == 0 || isnan(prodRateReal)) 343 // { 343 // { 344 // G4Exception("G4ReactionTableMessen 344 // G4Exception("G4ReactionTableMessenger", 345 // "WRONG_PRODUCTION_RATE 345 // "WRONG_PRODUCTION_RATE", 346 // FatalException, 346 // FatalException, 347 // ""); 347 // ""); 348 // } 348 // } 349 // 349 // 350 // double dimensionedProductionRate = p 350 // double dimensionedProductionRate = prodRateReal 351 // * (1e-3 * m3 / (mole * s)); 351 // * (1e-3 * m3 / (mole * s)); 352 // reactionData->SetProductionRate(dime 352 // reactionData->SetProductionRate(dimensionedProductionRate); 353 // } 353 // } 354 } 354 } 355 else if(rateconst_method == "Arr") 355 else if(rateconst_method == "Arr") 356 { 356 { 357 iss >> marker; // must be | 357 iss >> marker; // must be | 358 double A0 = 0; 358 double A0 = 0; 359 double E_R = 0; 359 double E_R = 0; 360 360 361 iss >> A0; 361 iss >> A0; 362 iss >> E_R; 362 iss >> E_R; 363 reactionData->SetArrehniusParameterizati 363 reactionData->SetArrehniusParameterization(A0, E_R); 364 } 364 } 365 else if(rateconst_method == "Pol") 365 else if(rateconst_method == "Pol") 366 { 366 { 367 iss >> marker; // must be | 367 iss >> marker; // must be | 368 std::vector<double> P = {0, 0, 0, 0, 0}; 368 std::vector<double> P = {0, 0, 0, 0, 0}; 369 369 370 size_t i = 0; 370 size_t i = 0; 371 while(i < 4) // could be changed to 5 on 371 while(i < 4) // could be changed to 5 only if marker is used as delimiter 372 { 372 { 373 double tmp; 373 double tmp; 374 iss >> tmp; 374 iss >> tmp; 375 P[i] = tmp; 375 P[i] = tmp; 376 // G4cout << newValue << G4endl; 376 // G4cout << newValue << G4endl; 377 // G4cout << tmp << G4endl; 377 // G4cout << tmp << G4endl; 378 // G4cout << P[i] << G4endl; 378 // G4cout << P[i] << G4endl; 379 ++i; 379 ++i; 380 }; 380 }; 381 reactionData->SetPolynomialParameterizat 381 reactionData->SetPolynomialParameterization(P); 382 } 382 } 383 else if(rateconst_method == "Scale") 383 else if(rateconst_method == "Scale") 384 { 384 { 385 iss >> marker; // must be | 385 iss >> marker; // must be | 386 double temp_K; 386 double temp_K; 387 iss >> temp_K; 387 iss >> temp_K; 388 double reactionRateCste; 388 double reactionRateCste; 389 iss >> reactionRateCste; 389 iss >> reactionRateCste; 390 double dimensionedReactionRate = reactio 390 double dimensionedReactionRate = reactionRateCste * (1e-3 * m3 / (mole * s)); 391 reactionData->SetObservedReactionRateCon 391 reactionData->SetObservedReactionRateConstant(dimensionedReactionRate); 392 reactionData->SetScaledParameterization( 392 reactionData->SetScaledParameterization(temp_K, dimensionedReactionRate); 393 } 393 } 394 394 395 // if(iss.eof() == false) 395 // if(iss.eof() == false) 396 // { 396 // { 397 // iss >> marker; 397 // iss >> marker; 398 // 398 // 399 // if(marker == "|") 399 // if(marker == "|") 400 // { 400 // { 401 // G4String productionRate ; 401 // G4String productionRate ; 402 // iss >> productionRate; 402 // iss >> productionRate; 403 // 403 // 404 //// G4cout << productionRate << G4endl 404 //// G4cout << productionRate << G4endl; 405 // 405 // 406 // double dimProductionRate = G4UIcomma 406 // double dimProductionRate = G4UIcommand::ConvertToDouble(productionRate)* (1e-3 * m3 / (mole * s)); 407 // 407 // 408 // G4cout << " DIM PROD RATE = " << rea 408 // G4cout << " DIM PROD RATE = " << reactionData->GetReactant1()->GetName() 409 // << " + " << reactionData->GetRe 409 // << " + " << reactionData->GetReactant2()->GetName() << " = " << dimProductionRate << G4endl; 410 // 410 // 411 // reactionData->SetProductionRate(dimP 411 // reactionData->SetProductionRate(dimProductionRate); 412 // } 412 // } 413 // } 413 // } 414 fpTable->SetReaction(reactionData); 414 fpTable->SetReaction(reactionData); 415 // G4cout << "Reaction " << species1 << " + 415 // G4cout << "Reaction " << species1 << " + " << species2 << " added" << G4endl; 416 } 416 } 417 } 417 } 418 418