60 : myA(0), myZ(0), theNucleons(250), currentNucleon(-1), theDensity(0),
61 nucleondistance(0.8*
fermi),excitationEnergy(0.),
62 places(250),
momentum(250), fermiM(250), testSums(250)
71 #if defined(NON_INTEGER_A_Z)
118 for (
G4int aNucleon=0; aNucleon <
myA; aNucleon++)
188 for (
int i=0; i<
myA; i++)
190 if (
theNucleons[i].GetPosition().mag2() > maxradius2 )
225 G4double factor=(1-std::sqrt(1-beta2))/beta2;
229 factor * (theBeta*
theNucleons[i].GetPosition()) * theBeta;
237 if (theBoost.
e() !=0 ) {
276 G4int protons=0,nucleons=0;
278 while (nucleons <
myA )
285 else if ( (nucleons-protons) < (
myA-
myZ) )
289 else G4cout <<
"G4Fancy3DNucleus::ChooseNucleons not efficient" <<
G4endl;
310 while ( (i <
myA) && (--interationsLeft>0))
317 G4RandFlat::shootArray(jr,prand);
322 aPos.
set((2*arand[jx]-1.), (2*arand[jy]-1.), (2*arand[--jr]-1.));
323 }
while (aPos.
mag2() > 1. );
329 std::vector<G4ThreeVector>::iterator iplace;
330 for( iplace=
places.begin(); iplace!=
places.end() && freeplace;++iplace)
332 delta = *iplace - aPos;
333 freeplace= delta.
mag2() > nd2;
343 G4double eFermi= std::sqrt(
sqr(pFermi) +
sqr(nucMass) ) - nucMass;
354 if (interationsLeft<=0) {
356 "Problem to place nucleons");
373 G4int loopCounterLeft = 10000;
374 for(
G4int ii=1; ii<4; ii++)
382 for(
G4int jj=0; jj < ii; jj++)
384 if( (
theNucleons[ii].GetPosition() -
theNucleons[jj].GetPosition()).mag2() <= nd2 ) {Continue =
true;
break;}
386 }
while( Continue && --loopCounterLeft > 0 );
388 if ( loopCounterLeft <= 0 ) {
390 "Unable to find a good position for the first alpha cluster");
392 loopCounterLeft = 10000;
393 for(
G4int ii=4; ii<8; ii++)
401 for(
G4int jj=0; jj < ii; jj++)
403 if( (
theNucleons[ii].GetPosition() -
theNucleons[jj].GetPosition()).mag2() <= nd2 ) {Continue =
true;
break;}
405 }
while( Continue && --loopCounterLeft > 0 );
407 if ( loopCounterLeft <= 0 ) {
409 "Unable to find a good position for the second alpha cluster");
411 loopCounterLeft = 10000;
412 for(
G4int ii=8; ii<12; ii++)
420 for(
G4int jj=0; jj < ii; jj++)
422 if( (
theNucleons[ii].GetPosition() -
theNucleons[jj].GetPosition()).mag2() <= nd2 ) {Continue =
true;
break;}
424 }
while( Continue && --loopCounterLeft > 0 );
426 if ( loopCounterLeft <= 0 ) {
428 "Unable to find a good position for the third alpha cluster");
454 for (
G4int ntry=0; ntry<1 ; ntry ++ )
456 for (i=0; i <
myA; i++ )
465 if ( eMax >
theNucleons[i].GetDefinition()->GetPDGMass() )
468 fermiM[i] = std::sqrt(pmax2);
469 while ( mom.
mag2() > pmax2 )
478 ed <<
"Nucleus Z A " <<
myZ <<
" " << myA <<
G4endl;
479 ed <<
"proton with eMax=" << eMax <<
G4endl;
480 G4Exception(
"G4Fancy3DNucleus::ChooseFermiMomenta(): difficulty finding proton momentum, set it to (0,0,0)",
499 for ( i=0; i<
myA ; i++ )
501 energy =
theNucleons[i].GetParticleType()->GetPDGMass()
521 if ( sum.
mag() <= PFermi )
533 for (
G4int aNucleon=0; aNucleon < myA-1; aNucleon++) {
534 delta = 2.*((
momentum[aNucleon]*testDir)*testDir);
536 testSums[aNucleon].Fill(delta, delta.
mag(), aNucleon);
543 while ( (sum-
testSums[--index].Vector).mag()>PFermi && index>0)
546 if ( sum.
mag() > (sum-
testSums[index].Vector).mag() ) {
548 sum-=testSums[index].Vector;
552 if ( (sum-
testSums[index].Vector).mag() <= PFermi )
556 for (
G4int aNucleon=0; aNucleon<=index; aNucleon++)
569 G4String text =
"G4Fancy3DNucleus.cc: Logic error in ReduceSum()";
583 if (
fermiM[++swapit] > PFermi )
break;
585 if (swapit == myA-1 )
return false;