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 // this code implementation is the intellectua 26 // this code implementation is the intellectual property of 27 // neutron_hp -- source file 27 // neutron_hp -- source file 28 // J.P. Wellisch, Nov-1996 28 // J.P. Wellisch, Nov-1996 29 // A prototype of the low energy neutron trans 29 // A prototype of the low energy neutron transport model. 30 // 30 // 31 // By copying, distributing or modifying the P 31 // By copying, distributing or modifying the Program (or any work 32 // based on the Program) you indicate your acc 32 // based on the Program) you indicate your acceptance of this statement, 33 // and all its terms. 33 // and all its terms. 34 // 34 // 35 // 35 // 36 // 070523 bug fix for G4FPE_DEBUG on by A. How 36 // 070523 bug fix for G4FPE_DEBUG on by A. Howard (and T. Koi) 37 // 081203 limit maximum trial for creating fin 37 // 081203 limit maximum trial for creating final states add protection for 1H isotope case by T. Koi 38 // 38 // 39 // P. Arce, June-2014 Conversion neutron_hp to 39 // P. Arce, June-2014 Conversion neutron_hp to particle_hp 40 // V. Ivanchenko, July-2023 Basic revision of << 41 // 40 // 42 #include "G4ParticleHPInelastic.hh" 41 #include "G4ParticleHPInelastic.hh" 43 << 42 #include "G4SystemOfUnits.hh" >> 43 #include "G4ParticleHPManager.hh" 44 #include "G4HadronicParameters.hh" 44 #include "G4HadronicParameters.hh" >> 45 #include "G4Threading.hh" >> 46 >> 47 G4ParticleHPInelastic::G4ParticleHPInelastic(G4ParticleDefinition* projectile, const char* name ) >> 48 :G4HadronicInteraction(name) >> 49 ,theInelastic(NULL) >> 50 ,numEle(0) >> 51 ,theProjectile(projectile) >> 52 { >> 53 G4String baseName; >> 54 if ( std::getenv("G4PARTICLEHPDATA") ) { >> 55 baseName = std::getenv( "G4PARTICLEHPDATA" ); >> 56 } >> 57 //const char* dataDirVariable; >> 58 G4String particleName; >> 59 if ( theProjectile == G4Neutron::Neutron() ) { >> 60 dataDirVariable = "G4NEUTRONHPDATA"; >> 61 } else if( theProjectile == G4Proton::Proton() ) { >> 62 dataDirVariable = "G4PROTONHPDATA"; >> 63 particleName = "Proton"; >> 64 } else if( theProjectile == G4Deuteron::Deuteron() ) { >> 65 dataDirVariable = "G4DEUTERONHPDATA"; >> 66 particleName = "Deuteron"; >> 67 } else if( theProjectile == G4Triton::Triton() ) { >> 68 dataDirVariable = "G4TRITONHPDATA"; >> 69 particleName = "Triton"; >> 70 } else if( theProjectile == G4He3::He3() ) { >> 71 dataDirVariable = "G4HE3HPDATA"; >> 72 particleName = "He3"; >> 73 } else if( theProjectile == G4Alpha::Alpha() ) { >> 74 dataDirVariable = "G4ALPHAHPDATA"; >> 75 particleName = "Alpha"; >> 76 } else { >> 77 G4String message("G4ParticleHPInelastic may only be called for neutron, proton, deuteron, triton, He3 or alpha, while it is called for " + theProjectile->GetParticleName()); >> 78 throw G4HadronicException(__FILE__, __LINE__,message.c_str()); >> 79 } >> 80 >> 81 SetMinEnergy( 0.0 ); >> 82 SetMaxEnergy( 20.*MeV ); >> 83 >> 84 //G4cout << " entering G4ParticleHPInelastic constructor"<<G4endl; >> 85 if ( !std::getenv("G4PARTICLEHPDATA") && !std::getenv(dataDirVariable) ) { >> 86 G4String message("Please setenv G4PARTICLEHPDATA (recommended) or, at least setenv " + >> 87 G4String(dataDirVariable) + " to point to the " + theProjectile->GetParticleName() + " cross-section files." ); >> 88 throw G4HadronicException(__FILE__, __LINE__,message.c_str()); >> 89 } >> 90 if ( std::getenv(dataDirVariable) ) { >> 91 dirName = std::getenv(dataDirVariable); >> 92 } else { >> 93 dirName = baseName + "/" + particleName; >> 94 } >> 95 >> 96 #ifdef G4VERBOSE >> 97 if ( G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) G4cout << dirName << G4endl; >> 98 #endif >> 99 >> 100 G4String tString = "/Inelastic"; >> 101 dirName = dirName + tString; >> 102 //numEle = G4Element::GetNumberOfElements(); >> 103 >> 104 #ifdef G4VERBOSE >> 105 if ( G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) >> 106 G4cout << "@@@ G4ParticleHPInelastic instantiated for particle " << theProjectile->GetParticleName() << " data directory variable is " << dataDirVariable << " pointing to " << dirName << G4endl; >> 107 #endif >> 108 >> 109 /* >> 110 theInelastic = new G4ParticleHPChannelList[numEle]; >> 111 for (G4int i=0; i<numEle; i++) >> 112 { >> 113 theInelastic[i].Init((*(G4Element::GetElementTable()))[i], dirName); >> 114 G4int itry = 0; >> 115 do >> 116 { >> 117 theInelastic[i].Register(&theNFS, "F01"); // has >> 118 theInelastic[i].Register(&theNXFS, "F02"); >> 119 theInelastic[i].Register(&the2NDFS, "F03"); >> 120 theInelastic[i].Register(&the2NFS, "F04"); // has, E Done >> 121 theInelastic[i].Register(&the3NFS, "F05"); // has, E Done >> 122 theInelastic[i].Register(&theNAFS, "F06"); >> 123 theInelastic[i].Register(&theN3AFS, "F07"); >> 124 theInelastic[i].Register(&the2NAFS, "F08"); >> 125 theInelastic[i].Register(&the3NAFS, "F09"); >> 126 theInelastic[i].Register(&theNPFS, "F10"); >> 127 theInelastic[i].Register(&theN2AFS, "F11"); >> 128 theInelastic[i].Register(&the2N2AFS, "F12"); >> 129 theInelastic[i].Register(&theNDFS, "F13"); >> 130 theInelastic[i].Register(&theNTFS, "F14"); >> 131 theInelastic[i].Register(&theNHe3FS, "F15"); >> 132 theInelastic[i].Register(&theND2AFS, "F16"); >> 133 theInelastic[i].Register(&theNT2AFS, "F17"); >> 134 theInelastic[i].Register(&the4NFS, "F18"); // has, E Done >> 135 theInelastic[i].Register(&the2NPFS, "F19"); >> 136 theInelastic[i].Register(&the3NPFS, "F20"); >> 137 theInelastic[i].Register(&theN2PFS, "F21"); >> 138 theInelastic[i].Register(&theNPAFS, "F22"); >> 139 theInelastic[i].Register(&thePFS, "F23"); >> 140 theInelastic[i].Register(&theDFS, "F24"); >> 141 theInelastic[i].Register(&theTFS, "F25"); >> 142 theInelastic[i].Register(&theHe3FS, "F26"); >> 143 theInelastic[i].Register(&theAFS, "F27"); >> 144 theInelastic[i].Register(&the2AFS, "F28"); >> 145 theInelastic[i].Register(&the3AFS, "F29"); >> 146 theInelastic[i].Register(&the2PFS, "F30"); >> 147 theInelastic[i].Register(&thePAFS, "F31"); >> 148 theInelastic[i].Register(&theD2AFS, "F32"); >> 149 theInelastic[i].Register(&theT2AFS, "F33"); >> 150 theInelastic[i].Register(&thePDFS, "F34"); >> 151 theInelastic[i].Register(&thePTFS, "F35"); >> 152 theInelastic[i].Register(&theDAFS, "F36"); >> 153 theInelastic[i].RestartRegistration(); >> 154 itry++; >> 155 } >> 156 //while(!theInelastic[i].HasDataInAnyFinalState()); >> 157 while( !theInelastic[i].HasDataInAnyFinalState() && itry < 6 ); >> 158 // 6 is corresponding to the value(5) of G4ParticleHPChannel. TK >> 159 >> 160 if ( itry == 6 ) >> 161 { >> 162 // No Final State at all. >> 163 G4bool exceptional = false; >> 164 if ( (*(G4Element::GetElementTable()))[i]->GetNumberOfIsotopes() == 1 ) >> 165 { >> 166 if ( (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetZ() == 1 && (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetN() == 1 ) exceptional = true; //1H >> 167 } >> 168 if ( !exceptional ) throw G4HadronicException(__FILE__, __LINE__, "Channel: Do not know what to do with this element"); >> 169 } >> 170 } >> 171 */ >> 172 /* >> 173 for (G4int i=0; i<numEle; i++) >> 174 { >> 175 theInelastic.push_back( new G4ParticleHPChannelList ); >> 176 (*theInelastic[i]).Init((*(G4Element::GetElementTable()))[i], dirName, theProjectile); >> 177 G4int itry = 0; >> 178 do >> 179 { >> 180 (*theInelastic[i]).Register(&theNFS, "F01"); // has >> 181 (*theInelastic[i]).Register(&theNXFS, "F02"); >> 182 (*theInelastic[i]).Register(&the2NDFS, "F03"); >> 183 (*theInelastic[i]).Register(&the2NFS, "F04"); // has, E Done >> 184 (*theInelastic[i]).Register(&the3NFS, "F05"); // has, E Done >> 185 (*theInelastic[i]).Register(&theNAFS, "F06"); >> 186 (*theInelastic[i]).Register(&theN3AFS, "F07"); >> 187 (*theInelastic[i]).Register(&the2NAFS, "F08"); >> 188 (*theInelastic[i]).Register(&the3NAFS, "F09"); >> 189 (*theInelastic[i]).Register(&theNPFS, "F10"); >> 190 (*theInelastic[i]).Register(&theN2AFS, "F11"); >> 191 (*theInelastic[i]).Register(&the2N2AFS, "F12"); >> 192 (*theInelastic[i]).Register(&theNDFS, "F13"); >> 193 (*theInelastic[i]).Register(&theNTFS, "F14"); >> 194 (*theInelastic[i]).Register(&theNHe3FS, "F15"); >> 195 (*theInelastic[i]).Register(&theND2AFS, "F16"); >> 196 (*theInelastic[i]).Register(&theNT2AFS, "F17"); >> 197 (*theInelastic[i]).Register(&the4NFS, "F18"); // has, E Done >> 198 (*theInelastic[i]).Register(&the2NPFS, "F19"); >> 199 (*theInelastic[i]).Register(&the3NPFS, "F20"); >> 200 (*theInelastic[i]).Register(&theN2PFS, "F21"); >> 201 (*theInelastic[i]).Register(&theNPAFS, "F22"); >> 202 (*theInelastic[i]).Register(&thePFS, "F23"); >> 203 (*theInelastic[i]).Register(&theDFS, "F24"); >> 204 (*theInelastic[i]).Register(&theTFS, "F25"); >> 205 (*theInelastic[i]).Register(&theHe3FS, "F26"); >> 206 (*theInelastic[i]).Register(&theAFS, "F27"); >> 207 (*theInelastic[i]).Register(&the2AFS, "F28"); >> 208 (*theInelastic[i]).Register(&the3AFS, "F29"); >> 209 (*theInelastic[i]).Register(&the2PFS, "F30"); >> 210 (*theInelastic[i]).Register(&thePAFS, "F31"); >> 211 (*theInelastic[i]).Register(&theD2AFS, "F32"); >> 212 (*theInelastic[i]).Register(&theT2AFS, "F33"); >> 213 (*theInelastic[i]).Register(&thePDFS, "F34"); >> 214 (*theInelastic[i]).Register(&thePTFS, "F35"); >> 215 (*theInelastic[i]).Register(&theDAFS, "F36"); >> 216 (*theInelastic[i]).RestartRegistration(); >> 217 itry++; >> 218 } >> 219 while( !(*theInelastic[i]).HasDataInAnyFinalState() && itry < 6 ); >> 220 // 6 is corresponding to the value(5) of G4ParticleHPChannel. TK >> 221 >> 222 if ( itry == 6 ) >> 223 { >> 224 // No Final State at all. >> 225 G4bool exceptional = false; >> 226 if ( (*(G4Element::GetElementTable()))[i]->GetNumberOfIsotopes() == 1 ) >> 227 { >> 228 if ( (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetZ() == 1 && (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetN() == 1 ) exceptional = true; //1H >> 229 } >> 230 if ( !exceptional ) { >> 231 G4cerr << " ELEMENT Z " << (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetZ() << " N " << (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetN() << G4endl; //1H >> 232 throw G4HadronicException(__FILE__, __LINE__, "Channel: Do not know what to do with this element"); >> 233 } >> 234 } >> 235 >> 236 } >> 237 */ >> 238 >> 239 } >> 240 >> 241 G4ParticleHPInelastic::~G4ParticleHPInelastic() >> 242 { >> 243 // delete [] theInelastic; >> 244 //Vector is shared, only master deletes >> 245 if ( !G4Threading::IsWorkerThread() ) { >> 246 if ( theInelastic != NULL ) { >> 247 for ( std::vector<G4ParticleHPChannelList*>::iterator >> 248 it = theInelastic->begin() ; it != theInelastic->end() ; it++ ) { >> 249 delete *it; >> 250 } >> 251 theInelastic->clear(); >> 252 } >> 253 } >> 254 } >> 255 >> 256 #include "G4ParticleHPThermalBoost.hh" >> 257 >> 258 G4HadFinalState * G4ParticleHPInelastic::ApplyYourself(const G4HadProjectile& aTrack, G4Nucleus& aNucleus ) >> 259 { >> 260 G4ParticleHPManager::GetInstance()->OpenReactionWhiteBoard(); >> 261 const G4Material * theMaterial = aTrack.GetMaterial(); >> 262 G4int n = theMaterial->GetNumberOfElements(); >> 263 G4int index = theMaterial->GetElement(0)->GetIndex(); >> 264 G4int it=0; >> 265 if(n!=1) >> 266 { >> 267 G4double* xSec = new G4double[n]; >> 268 G4double sum=0; >> 269 G4int i; >> 270 const G4double * NumAtomsPerVolume = theMaterial->GetVecNbOfAtomsPerVolume(); >> 271 G4double rWeight; >> 272 G4ParticleHPThermalBoost aThermalE; >> 273 for (i=0; i<n; i++) >> 274 { >> 275 index = theMaterial->GetElement(i)->GetIndex(); >> 276 rWeight = NumAtomsPerVolume[i]; >> 277 if ( aTrack.GetDefinition() == G4Neutron::Neutron() ) { >> 278 xSec[i] = ((*theInelastic)[index])->GetXsec(aThermalE.GetThermalEnergy(aTrack, >> 279 theMaterial->GetElement(i), >> 280 theMaterial->GetTemperature())); >> 281 } else { >> 282 xSec[i] = ((*theInelastic)[index])->GetXsec(aTrack.GetKineticEnergy()); >> 283 } >> 284 xSec[i] *= rWeight; >> 285 sum+=xSec[i]; >> 286 #ifdef G4PHPDEBUG >> 287 #ifdef G4VERBOSE >> 288 if( std::getenv("G4ParticleHPDebug") && G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) >> 289 G4cout << " G4ParticleHPInelastic XSEC ELEM " << i << " = " << xSec[i] << G4endl; >> 290 #endif >> 291 #endif >> 292 >> 293 } >> 294 G4double random = G4UniformRand(); >> 295 G4double running = 0; >> 296 for (i=0; i<n; i++) >> 297 { >> 298 running += xSec[i]; >> 299 index = theMaterial->GetElement(i)->GetIndex(); >> 300 it = i; >> 301 //if(random<=running/sum) break; >> 302 if( sum == 0 || random<=running/sum) break; >> 303 } >> 304 delete [] xSec; >> 305 } >> 306 >> 307 #ifdef G4PHPDEBUG >> 308 #ifdef G4VERBOSE >> 309 if( std::getenv("G4ParticleHPDebug") && G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) >> 310 G4cout << " G4ParticleHPInelastic SELECTED ELEM " << it << " = " << theMaterial->GetElement(it)->GetName() << " FROM MATERIAL " << theMaterial->GetName() << G4endl; >> 311 #endif >> 312 #endif >> 313 //return theInelastic[index].ApplyYourself(theMaterial->GetElement(it), aTrack); >> 314 G4HadFinalState* result = ((*theInelastic)[index])->ApplyYourself(theMaterial->GetElement(it), aTrack); >> 315 // >> 316 aNucleus.SetParameters(G4ParticleHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargA(),G4ParticleHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargZ()); >> 317 const G4Element* target_element = (*G4Element::GetElementTable())[index]; >> 318 const G4Isotope* target_isotope=NULL; >> 319 G4int iele = target_element->GetNumberOfIsotopes(); >> 320 for ( G4int j = 0 ; j != iele ; j++ ) { >> 321 target_isotope=target_element->GetIsotope( j ); >> 322 if ( target_isotope->GetN() == G4ParticleHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargA() ) break; >> 323 } >> 324 //G4cout << "Target Material of this reaction is " << theMaterial->GetName() << G4endl; >> 325 //G4cout << "Target Element of this reaction is " << target_element->GetName() << G4endl; >> 326 //G4cout << "Target Isotope of this reaction is " << target_isotope->GetName() << G4endl; >> 327 aNucleus.SetIsotope( target_isotope ); >> 328 >> 329 G4ParticleHPManager::GetInstance()->CloseReactionWhiteBoard(); >> 330 >> 331 //GDEB >> 332 if( std::getenv("G4PHPTEST") ) { >> 333 G4HadSecondary* seco = result->GetSecondary(0); >> 334 if(seco) { >> 335 G4ThreeVector secoMom = seco->GetParticle()->GetMomentum(); >> 336 #ifdef G4VERBOSE >> 337 if ( G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) >> 338 G4cout << " G4ParticleHPinelastic COS THETA " << std::cos(secoMom.theta()) <<" " << secoMom << G4endl; >> 339 #endif >> 340 } >> 341 } >> 342 >> 343 return result; >> 344 } >> 345 >> 346 const std::pair<G4double, G4double> G4ParticleHPInelastic::GetFatalEnergyCheckLevels() const >> 347 { >> 348 // max energy non-conservation is mass of heavy nucleus >> 349 // This should be same to the hadron default value >> 350 return std::pair<G4double, G4double>(10*perCent,DBL_MAX); >> 351 } >> 352 >> 353 /* >> 354 void G4ParticleHPInelastic::addChannelForNewElement() >> 355 { >> 356 for ( G4int i = numEle ; i < (G4int)G4Element::GetNumberOfElements() ; i++ ) >> 357 { >> 358 G4cout << "G4ParticleHPInelastic Prepairing Data for the new element of " << (*(G4Element::GetElementTable()))[i]->GetName() << G4endl; >> 359 >> 360 theInelastic.push_back( new G4ParticleHPChannelList ); >> 361 (*theInelastic[i]).Init((*(G4Element::GetElementTable()))[i], dirName, theProjectile); >> 362 G4int itry = 0; >> 363 do >> 364 { >> 365 (*theInelastic[i]).Register(&theNFS, "F01"); // has >> 366 (*theInelastic[i]).Register(&theNXFS, "F02"); >> 367 (*theInelastic[i]).Register(&the2NDFS, "F03"); >> 368 (*theInelastic[i]).Register(&the2NFS, "F04"); // has, E Done >> 369 (*theInelastic[i]).Register(&the3NFS, "F05"); // has, E Done >> 370 (*theInelastic[i]).Register(&theNAFS, "F06"); >> 371 (*theInelastic[i]).Register(&theN3AFS, "F07"); >> 372 (*theInelastic[i]).Register(&the2NAFS, "F08"); >> 373 (*theInelastic[i]).Register(&the3NAFS, "F09"); >> 374 (*theInelastic[i]).Register(&theNPFS, "F10"); >> 375 (*theInelastic[i]).Register(&theN2AFS, "F11"); >> 376 (*theInelastic[i]).Register(&the2N2AFS, "F12"); >> 377 (*theInelastic[i]).Register(&theNDFS, "F13"); >> 378 (*theInelastic[i]).Register(&theNTFS, "F14"); >> 379 (*theInelastic[i]).Register(&theNHe3FS, "F15"); >> 380 (*theInelastic[i]).Register(&theND2AFS, "F16"); >> 381 (*theInelastic[i]).Register(&theNT2AFS, "F17"); >> 382 (*theInelastic[i]).Register(&the4NFS, "F18"); // has, E Done >> 383 (*theInelastic[i]).Register(&the2NPFS, "F19"); >> 384 (*theInelastic[i]).Register(&the3NPFS, "F20"); >> 385 (*theInelastic[i]).Register(&theN2PFS, "F21"); >> 386 (*theInelastic[i]).Register(&theNPAFS, "F22"); >> 387 (*theInelastic[i]).Register(&thePFS, "F23"); >> 388 (*theInelastic[i]).Register(&theDFS, "F24"); >> 389 (*theInelastic[i]).Register(&theTFS, "F25"); >> 390 (*theInelastic[i]).Register(&theHe3FS, "F26"); >> 391 (*theInelastic[i]).Register(&theAFS, "F27"); >> 392 (*theInelastic[i]).Register(&the2AFS, "F28"); >> 393 (*theInelastic[i]).Register(&the3AFS, "F29"); >> 394 (*theInelastic[i]).Register(&the2PFS, "F30"); >> 395 (*theInelastic[i]).Register(&thePAFS, "F31"); >> 396 (*theInelastic[i]).Register(&theD2AFS, "F32"); >> 397 (*theInelastic[i]).Register(&theT2AFS, "F33"); >> 398 (*theInelastic[i]).Register(&thePDFS, "F34"); >> 399 (*theInelastic[i]).Register(&thePTFS, "F35"); >> 400 (*theInelastic[i]).Register(&theDAFS, "F36"); >> 401 (*theInelastic[i]).RestartRegistration(); >> 402 itry++; >> 403 } >> 404 while( !(*theInelastic[i]).HasDataInAnyFinalState() && itry < 6 ); >> 405 // 6 is corresponding to the value(5) of G4ParticleHPChannel. TK >> 406 >> 407 if ( itry == 6 ) >> 408 { >> 409 // No Final State at all. >> 410 G4bool exceptional = false; >> 411 if ( (*(G4Element::GetElementTable()))[i]->GetNumberOfIsotopes() == 1 ) >> 412 { >> 413 if ( (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetZ() == 1 && (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetN() == 1 ) exceptional = true; //1H >> 414 } >> 415 if ( !exceptional ) throw G4HadronicException(__FILE__, __LINE__, "Channel: Do not know what to do with this element"); >> 416 } >> 417 } >> 418 >> 419 numEle = (G4int)G4Element::GetNumberOfElements(); >> 420 } >> 421 */ >> 422 G4int G4ParticleHPInelastic::GetVerboseLevel() const >> 423 { >> 424 return G4ParticleHPManager::GetInstance()->GetVerboseLevel(); >> 425 } >> 426 void G4ParticleHPInelastic::SetVerboseLevel( G4int newValue ) >> 427 { >> 428 G4ParticleHPManager::GetInstance()->SetVerboseLevel(newValue); >> 429 } >> 430 45 #include "G4ParticleHP2AInelasticFS.hh" 431 #include "G4ParticleHP2AInelasticFS.hh" 46 #include "G4ParticleHP2N2AInelasticFS.hh" 432 #include "G4ParticleHP2N2AInelasticFS.hh" 47 #include "G4ParticleHP2NAInelasticFS.hh" 433 #include "G4ParticleHP2NAInelasticFS.hh" 48 #include "G4ParticleHP2NDInelasticFS.hh" 434 #include "G4ParticleHP2NDInelasticFS.hh" 49 #include "G4ParticleHP2NInelasticFS.hh" 435 #include "G4ParticleHP2NInelasticFS.hh" 50 #include "G4ParticleHP2NPInelasticFS.hh" 436 #include "G4ParticleHP2NPInelasticFS.hh" 51 #include "G4ParticleHP2PInelasticFS.hh" 437 #include "G4ParticleHP2PInelasticFS.hh" 52 #include "G4ParticleHP3AInelasticFS.hh" 438 #include "G4ParticleHP3AInelasticFS.hh" 53 #include "G4ParticleHP3NAInelasticFS.hh" 439 #include "G4ParticleHP3NAInelasticFS.hh" 54 #include "G4ParticleHP3NInelasticFS.hh" 440 #include "G4ParticleHP3NInelasticFS.hh" 55 #include "G4ParticleHP3NPInelasticFS.hh" 441 #include "G4ParticleHP3NPInelasticFS.hh" 56 #include "G4ParticleHP4NInelasticFS.hh" 442 #include "G4ParticleHP4NInelasticFS.hh" 57 #include "G4ParticleHPAInelasticFS.hh" 443 #include "G4ParticleHPAInelasticFS.hh" 58 #include "G4ParticleHPD2AInelasticFS.hh" 444 #include "G4ParticleHPD2AInelasticFS.hh" 59 #include "G4ParticleHPDAInelasticFS.hh" 445 #include "G4ParticleHPDAInelasticFS.hh" 60 #include "G4ParticleHPDInelasticFS.hh" 446 #include "G4ParticleHPDInelasticFS.hh" 61 #include "G4ParticleHPHe3InelasticFS.hh" 447 #include "G4ParticleHPHe3InelasticFS.hh" 62 #include "G4ParticleHPManager.hh" << 63 #include "G4ParticleHPN2AInelasticFS.hh" 448 #include "G4ParticleHPN2AInelasticFS.hh" 64 #include "G4ParticleHPN2PInelasticFS.hh" 449 #include "G4ParticleHPN2PInelasticFS.hh" 65 #include "G4ParticleHPN3AInelasticFS.hh" 450 #include "G4ParticleHPN3AInelasticFS.hh" 66 #include "G4ParticleHPNAInelasticFS.hh" 451 #include "G4ParticleHPNAInelasticFS.hh" 67 #include "G4ParticleHPND2AInelasticFS.hh" 452 #include "G4ParticleHPND2AInelasticFS.hh" 68 #include "G4ParticleHPNDInelasticFS.hh" 453 #include "G4ParticleHPNDInelasticFS.hh" 69 #include "G4ParticleHPNHe3InelasticFS.hh" 454 #include "G4ParticleHPNHe3InelasticFS.hh" 70 #include "G4ParticleHPNInelasticFS.hh" 455 #include "G4ParticleHPNInelasticFS.hh" 71 #include "G4ParticleHPNPAInelasticFS.hh" 456 #include "G4ParticleHPNPAInelasticFS.hh" 72 #include "G4ParticleHPNPInelasticFS.hh" 457 #include "G4ParticleHPNPInelasticFS.hh" 73 #include "G4ParticleHPNT2AInelasticFS.hh" 458 #include "G4ParticleHPNT2AInelasticFS.hh" 74 #include "G4ParticleHPNTInelasticFS.hh" 459 #include "G4ParticleHPNTInelasticFS.hh" 75 #include "G4ParticleHPNXInelasticFS.hh" 460 #include "G4ParticleHPNXInelasticFS.hh" 76 #include "G4ParticleHPPAInelasticFS.hh" 461 #include "G4ParticleHPPAInelasticFS.hh" 77 #include "G4ParticleHPPDInelasticFS.hh" 462 #include "G4ParticleHPPDInelasticFS.hh" 78 #include "G4ParticleHPPInelasticFS.hh" 463 #include "G4ParticleHPPInelasticFS.hh" 79 #include "G4ParticleHPPTInelasticFS.hh" 464 #include "G4ParticleHPPTInelasticFS.hh" 80 #include "G4ParticleHPT2AInelasticFS.hh" 465 #include "G4ParticleHPT2AInelasticFS.hh" 81 #include "G4ParticleHPTInelasticFS.hh" 466 #include "G4ParticleHPTInelasticFS.hh" 82 #include "G4ParticleHPThermalBoost.hh" << 83 #include "G4SystemOfUnits.hh" << 84 #include "G4AutoLock.hh" << 85 467 86 G4bool G4ParticleHPInelastic::fLock[] = {true, << 468 void G4ParticleHPInelastic::BuildPhysicsTable(const G4ParticleDefinition& projectile) { 87 std::vector<G4ParticleHPChannelList*>* << 88 G4ParticleHPInelastic::theInelastic[] = {nullp << 89 << 90 namespace << 91 { << 92 G4Mutex theHPInelastic = G4MUTEX_INITIALIZER << 93 } << 94 469 95 G4ParticleHPInelastic::G4ParticleHPInelastic(G << 470 G4ParticleHPManager* hpmanager = G4ParticleHPManager::GetInstance(); 96 : G4HadronicInteraction(name), theProjectile << 97 { << 98 fManager = G4ParticleHPManager::GetInstance( << 99 dirName = fManager->GetParticleHPPath(thePro << 100 indexP = fManager->GetPHPIndex(theProjectile << 101 << 102 #ifdef G4VERBOSE << 103 if (fManager->GetVerboseLevel() > 1) << 104 G4cout << "@@@ G4ParticleHPInelastic insta << 105 << p->GetParticleName() << " indexP << 106 << "/n data directory " << dirName << << 107 #endif << 108 } << 109 << 110 G4ParticleHPInelastic::~G4ParticleHPInelastic( << 111 { << 112 // Vector is shared, only one delete << 113 if (isFirst) { << 114 ClearData(); << 115 } << 116 } << 117 471 118 void G4ParticleHPInelastic::ClearData() << 472 theInelastic = hpmanager->GetInelasticFinalStates( &projectile ); 119 { << 120 if (theInelastic[indexP] != nullptr) { << 121 for (auto const& p : *(theInelastic[indexP << 122 delete p; << 123 } << 124 delete theInelastic[indexP]; << 125 theInelastic[indexP] = nullptr; << 126 } << 127 } << 128 473 129 G4HadFinalState* G4ParticleHPInelastic::ApplyY << 474 if ( G4Threading::IsMasterThread() ) { 130 << 131 { << 132 G4ParticleHPManager::GetInstance()->OpenReac << 133 const G4Material* theMaterial = aTrack.GetMa << 134 auto n = (G4int)theMaterial->GetNumberOfElem << 135 auto elm = theMaterial->GetElement(0); << 136 std::size_t index = elm->GetIndex(); << 137 G4int it = 0; << 138 /* << 139 G4cout << "G4ParticleHPInelastic::ApplyYours << 140 << " indexP=" << indexP << " " << 141 << aTrack.GetDefinition()->GetParticl << 142 */ << 143 if (n != 1) { << 144 auto xSec = new G4double[n]; << 145 G4double sum = 0; << 146 G4int i; << 147 const G4double* NumAtomsPerVolume = theMat << 148 G4double rWeight; << 149 G4double xs; << 150 G4ParticleHPThermalBoost aThermalE; << 151 for (i = 0; i < n; ++i) { << 152 elm = theMaterial->GetElement(i); << 153 index = elm->GetIndex(); << 154 /* << 155 G4cout << "i=" << i << " index=" << ind << 156 << " " << (*(theInelastic[indexP]))[in << 157 */ << 158 rWeight = NumAtomsPerVolume[i]; << 159 if (aTrack.GetDefinition() == G4Neutron: << 160 xs = (*(theInelastic[indexP]))[index]- << 161 theMaterial->GetTemperature << 162 } << 163 else { << 164 xs = (*(theInelastic[indexP]))[index]- << 165 } << 166 xs *= rWeight; << 167 sum += xs; << 168 xSec[i] = sum; << 169 #ifdef G4VERBOSE << 170 if (fManager->GetDEBUG()) << 171 G4cout << " G4ParticleHPInelastic XSEC << 172 #endif << 173 } << 174 sum *= G4UniformRand(); << 175 for (it = 0; it < n; ++it) { << 176 elm = theMaterial->GetElement(it); << 177 index = elm->GetIndex(); << 178 if (sum <= xSec[it]) break; << 179 } << 180 delete[] xSec; << 181 } << 182 475 183 #ifdef G4VERBOSE << 476 if ( theInelastic == NULL ) theInelastic = new std::vector<G4ParticleHPChannelList*>; 184 if (fManager->GetDEBUG()) << 185 G4cout << " G4ParticleHPInelastic: Elem it << 186 << elm->GetName() << " index=" << i << 187 << " from material " << theMaterial->GetN << 188 << G4endl; << 189 #endif << 190 477 191 G4HadFinalState* result = << 478 if ( numEle == (G4int)G4Element::GetNumberOfElements() ) return; 192 (*(theInelastic[indexP]))[index]->ApplyYou << 193 479 194 aNucleus.SetParameters(fManager->GetReaction << 480 if ( theInelastic->size() == G4Element::GetNumberOfElements() ) { 195 fManager->GetReaction << 481 numEle = G4Element::GetNumberOfElements(); 196 << 482 return; 197 const G4Element* target_element = (*G4Elemen << 483 } 198 const G4Isotope* target_isotope = nullptr; << 484 /* 199 auto iele = (G4int)target_element->GetNumber << 485 const char* dataDirVariable; 200 for (G4int j = 0; j != iele; ++j) { << 486 if( &projectile == G4Neutron::Neutron() ) { 201 target_isotope = target_element->GetIsotop << 487 dataDirVariable = "G4NEUTRONHPDATA"; 202 if (target_isotope->GetN() == fManager->Ge << 488 }else if( &projectile == G4Proton::Proton() ) { 203 break; << 489 dataDirVariable = "G4PROTONHPDATA"; 204 } << 490 }else if( &projectile == G4Deuteron::Deuteron() ) { 205 aNucleus.SetIsotope(target_isotope); << 491 dataDirVariable = "G4DEUTERONHPDATA"; 206 << 492 }else if( &projectile == G4Triton::Triton() ) { 207 G4ParticleHPManager::GetInstance()->CloseRea << 493 dataDirVariable = "G4TRITONHPDATA"; 208 << 494 }else if( &projectile == G4He3::He3() ) { 209 return result; << 495 dataDirVariable = "G4HE3HPDATA"; 210 } << 496 }else if( &projectile == G4Alpha::Alpha() ) { 211 << 497 dataDirVariable = "G4ALPHAHPDATA"; 212 const std::pair<G4double, G4double> G4Particle << 498 } else { 213 { << 499 G4String message("G4ParticleHPInelastic may only be called for neutron, proton, deuteron, triton, He3 or alpha, while it is called for " + projectile.GetParticleName()); 214 // max energy non-conservation is mass of he << 500 throw G4HadronicException(__FILE__, __LINE__,message.c_str()); 215 return std::pair<G4double, G4double>(10.0 * << 501 } 216 } << 502 if(!std::getenv(dataDirVariable)){ 217 << 503 G4String message("Please set the environment variable " + G4String(dataDirVariable) + " to point to the " + projectile.GetParticleName() + " cross-section files."); 218 void G4ParticleHPInelastic::BuildPhysicsTable( << 504 throw G4HadronicException(__FILE__, __LINE__,message.c_str()); 219 { << 505 } 220 if (fLock[indexP]) { << 506 dirName = std::getenv(dataDirVariable); 221 G4AutoLock l(&theHPInelastic); << 507 G4cout << dirName << G4endl; 222 if (fLock[indexP]) { << 223 isFirst = true; << 224 fLock[indexP] = false; << 225 } << 226 l.unlock(); << 227 } << 228 << 229 G4int nelm = (G4int)G4Element::GetNumberOfEl << 230 G4int n0 = numEle; << 231 numEle = nelm; << 232 if (!isFirst || nelm == n0) { return; } << 233 << 234 // extra elements should be initialized << 235 G4AutoLock l(&theHPInelastic); << 236 << 237 if (nullptr == theInelastic[indexP]) { << 238 theInelastic[indexP] = new std::vector<G4P << 239 } << 240 508 241 if (fManager->GetVerboseLevel() > 0 && isFir << 509 G4String tString = "/Inelastic"; 242 fManager->DumpSetting(); << 510 dirName = dirName + tString; 243 G4cout << "@@@ G4ParticleHPInelastic insta << 244 << theProjectile->GetParticleName() << "/ << 245 << dirName << G4endl; << 246 } << 247 511 248 auto table = G4Element::GetElementTable(); << 512 */ 249 for (G4int i = n0; i < nelm; ++i) { << 513 #ifdef G4VERBOSE 250 auto clist = new G4ParticleHPChannelList(3 << 514 if ( G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) { 251 theInelastic[indexP]->push_back(clist); << 515 hpmanager->DumpSetting(); 252 clist->Init((*table)[i], dirName, theProje << 516 G4cout << "@@@ G4ParticleHPInelastic instantiated for particle " << projectile.GetParticleName() << " data directory variable is " << dataDirVariable << " pointing to " << dirName << G4endl; 253 clist->Register(new G4ParticleHPNInelastic << 517 } 254 clist->Register(new G4ParticleHPNXInelasti << 518 #endif 255 clist->Register(new G4ParticleHP2NDInelast << 519 for (G4int i = numEle ; i < (G4int)G4Element::GetNumberOfElements(); i++) 256 clist->Register(new G4ParticleHP2NInelasti << 520 { 257 clist->Register(new G4ParticleHP3NInelasti << 521 theInelastic->push_back( new G4ParticleHPChannelList ); 258 clist->Register(new G4ParticleHPNAInelasti << 522 ((*theInelastic)[i])->Init((*(G4Element::GetElementTable()))[i], dirName, const_cast<G4ParticleDefinition*>(&projectile)); 259 clist->Register(new G4ParticleHPN3AInelast << 523 G4int itry = 0; 260 clist->Register(new G4ParticleHP2NAInelast << 524 do 261 clist->Register(new G4ParticleHP3NAInelast << 525 { 262 clist->Register(new G4ParticleHPNPInelasti << 526 ((*theInelastic)[i])->Register( new G4ParticleHPNInelasticFS , "F01"); // has 263 clist->Register(new G4ParticleHPN2AInelast << 527 ((*theInelastic)[i])->Register( new G4ParticleHPNXInelasticFS , "F02"); 264 clist->Register(new G4ParticleHP2N2AInelas << 528 ((*theInelastic)[i])->Register( new G4ParticleHP2NDInelasticFS , "F03"); 265 clist->Register(new G4ParticleHPNDInelasti << 529 ((*theInelastic)[i])->Register( new G4ParticleHP2NInelasticFS , "F04"); // has, E Done 266 clist->Register(new G4ParticleHPNTInelasti << 530 ((*theInelastic)[i])->Register( new G4ParticleHP3NInelasticFS , "F05"); // has, E Done 267 clist->Register(new G4ParticleHPNHe3Inelas << 531 ((*theInelastic)[i])->Register( new G4ParticleHPNAInelasticFS , "F06"); 268 clist->Register(new G4ParticleHPND2AInelas << 532 ((*theInelastic)[i])->Register( new G4ParticleHPN3AInelasticFS , "F07"); 269 clist->Register(new G4ParticleHPNT2AInelas << 533 ((*theInelastic)[i])->Register( new G4ParticleHP2NAInelasticFS , "F08"); 270 clist->Register(new G4ParticleHP4NInelasti << 534 ((*theInelastic)[i])->Register( new G4ParticleHP3NAInelasticFS , "F09"); 271 clist->Register(new G4ParticleHP2NPInelast << 535 ((*theInelastic)[i])->Register( new G4ParticleHPNPInelasticFS , "F10"); 272 clist->Register(new G4ParticleHP3NPInelast << 536 ((*theInelastic)[i])->Register( new G4ParticleHPN2AInelasticFS , "F11"); 273 clist->Register(new G4ParticleHPN2PInelast << 537 ((*theInelastic)[i])->Register( new G4ParticleHP2N2AInelasticFS , "F12"); 274 clist->Register(new G4ParticleHPNPAInelast << 538 ((*theInelastic)[i])->Register( new G4ParticleHPNDInelasticFS , "F13"); 275 clist->Register(new G4ParticleHPPInelastic << 539 ((*theInelastic)[i])->Register( new G4ParticleHPNTInelasticFS , "F14"); 276 clist->Register(new G4ParticleHPDInelastic << 540 ((*theInelastic)[i])->Register( new G4ParticleHPNHe3InelasticFS , "F15"); 277 clist->Register(new G4ParticleHPTInelastic << 541 ((*theInelastic)[i])->Register( new G4ParticleHPND2AInelasticFS , "F16"); 278 clist->Register(new G4ParticleHPHe3Inelast << 542 ((*theInelastic)[i])->Register( new G4ParticleHPNT2AInelasticFS , "F17"); 279 clist->Register(new G4ParticleHPAInelastic << 543 ((*theInelastic)[i])->Register( new G4ParticleHP4NInelasticFS , "F18"); // has, E Done 280 clist->Register(new G4ParticleHP2AInelasti << 544 ((*theInelastic)[i])->Register( new G4ParticleHP2NPInelasticFS , "F19"); 281 clist->Register(new G4ParticleHP3AInelasti << 545 ((*theInelastic)[i])->Register( new G4ParticleHP3NPInelasticFS , "F20"); 282 clist->Register(new G4ParticleHP2PInelasti << 546 ((*theInelastic)[i])->Register( new G4ParticleHPN2PInelasticFS , "F21"); 283 clist->Register(new G4ParticleHPPAInelasti << 547 ((*theInelastic)[i])->Register( new G4ParticleHPNPAInelasticFS , "F22"); 284 clist->Register(new G4ParticleHPD2AInelast << 548 ((*theInelastic)[i])->Register( new G4ParticleHPPInelasticFS , "F23"); 285 clist->Register(new G4ParticleHPT2AInelast << 549 ((*theInelastic)[i])->Register( new G4ParticleHPDInelasticFS , "F24"); 286 clist->Register(new G4ParticleHPPDInelasti << 550 ((*theInelastic)[i])->Register( new G4ParticleHPTInelasticFS , "F25"); 287 clist->Register(new G4ParticleHPPTInelasti << 551 ((*theInelastic)[i])->Register( new G4ParticleHPHe3InelasticFS , "F26"); 288 clist->Register(new G4ParticleHPDAInelasti << 552 ((*theInelastic)[i])->Register( new G4ParticleHPAInelasticFS , "F27"); 289 #ifdef G4VERBOSE << 553 ((*theInelastic)[i])->Register( new G4ParticleHP2AInelasticFS , "F28"); 290 if (fManager->GetVerboseLevel() > 1) { << 554 ((*theInelastic)[i])->Register( new G4ParticleHP3AInelasticFS , "F29"); 291 G4cout << "ParticleHP::Inelastic for " << 555 ((*theInelastic)[i])->Register( new G4ParticleHP2PInelasticFS , "F30"); 292 << theProjectile->GetParticleName() << << 556 ((*theInelastic)[i])->Register( new G4ParticleHPPAInelasticFS , "F31"); 293 << (*table)[i]->GetName() << G4endl; << 557 ((*theInelastic)[i])->Register( new G4ParticleHPD2AInelasticFS , "F32"); 294 } << 558 ((*theInelastic)[i])->Register( new G4ParticleHPT2AInelasticFS , "F33"); 295 #endif << 559 ((*theInelastic)[i])->Register( new G4ParticleHPPDInelasticFS , "F34"); 296 } << 560 ((*theInelastic)[i])->Register( new G4ParticleHPPTInelasticFS , "F35"); 297 fManager->RegisterInelasticFinalStates( &pro << 561 ((*theInelastic)[i])->Register( new G4ParticleHPDAInelasticFS , "F36"); 298 l.unlock(); << 562 ((*theInelastic)[i])->RestartRegistration(); >> 563 itry++; >> 564 } >> 565 while( !((*theInelastic)[i])->HasDataInAnyFinalState() && itry < 6 ); // Loop checking, 11.05.2015, T. Koi >> 566 // 6 is corresponding to the value(5) of G4ParticleHPChannel. TK >> 567 >> 568 if ( itry == 6 ) { >> 569 // No Final State at all. >> 570 /* >> 571 G4bool exceptional = false; >> 572 if ( (*(G4Element::GetElementTable()))[i]->GetNumberOfIsotopes() == 1 ) >> 573 { >> 574 if ( (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetZ() == 1 && (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetN() == 1 ) exceptional = true; //1H >> 575 } >> 576 if ( !exceptional ) { >> 577 G4cerr << " ELEMENT Z " << (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetZ() << " N " << (*(G4Element::GetElementTable()))[i]->GetIsotope( 0 )->GetN() << G4endl; //1H >> 578 throw G4HadronicException(__FILE__, __LINE__, "Channel: Do not know what to do with this element"); >> 579 } >> 580 */ >> 581 #ifdef G4VERBOSE >> 582 if ( G4ParticleHPManager::GetInstance()->GetVerboseLevel() > 1 && >> 583 G4HadronicParameters::Instance()->GetVerboseLevel() > 0 ) { >> 584 G4cout << "ParticleHP::Inelastic for " << projectile.GetParticleName() << ". Do not know what to do with element of \"" << (*(G4Element::GetElementTable()))[i]->GetName() << "\"." << G4endl; >> 585 G4cout << "The components of the element are" << G4endl; >> 586 //G4cout << "TKDB dataDirVariable = " << dataDirVariable << G4endl; >> 587 for ( G4int ii = 0 ; ii < (G4int)( (*(G4Element::GetElementTable()))[i]->GetNumberOfIsotopes() ) ; ii++ ) { >> 588 G4cout << " Z: " << (*(G4Element::GetElementTable()))[i]->GetIsotope( ii )->GetZ() << ", A: " << (*(G4Element::GetElementTable()))[i]->GetIsotope( ii )->GetN() << G4endl; >> 589 } >> 590 G4cout << "No possible final state data of the element is found in " << dataDirVariable << "." << G4endl; >> 591 } >> 592 #endif >> 593 } >> 594 } >> 595 hpmanager->RegisterInelasticFinalStates( &projectile , theInelastic ); >> 596 } >> 597 numEle = G4Element::GetNumberOfElements(); 299 } 598 } 300 599 301 void G4ParticleHPInelastic::ModelDescription(s 600 void G4ParticleHPInelastic::ModelDescription(std::ostream& outFile) const 302 { 601 { 303 outFile << "High Precision (HP) model for in << 602 outFile << "Extension of High Precision model for inelastic reaction of proton, deuteron, triton, He3 and alpha below 20MeV\n"; 304 << theProjectile->GetParticleName() << " b << 305 } 603 } >> 604 306 605