Problem 1629

Summary: Bug in G4PhaseSpaceDecayChannel
Product: Geant4 Reporter: Witold.Pokorski
Component: particlesAssignee: kurasige
Status: RESOLVED FIXED    
Severity: major CC: Alberto.Ribon, Federico.Carminati, kurasige
Priority: P5    
Version: 10.0   
Hardware: All   
OS: All   

Description Witold.Pokorski 2014-05-19 16:07:26 CEST
When decaying anti_lambda(1800) Geant4 always outputs the following message:

G4PhaseSpaceDecayChannel::TwoBodyDecayIt sum of daughter mass is larger than parent mass
parent :anti_lambda(1800)  1.8
daughter 1 :anti_neutron  0.939565
daughter 2:k_star0  0.89594

and does not produce any decay products.

The reason for this is that there is an inconsistency between what is implemented in G4VDecayChannel and what is in G4PhaseSpaceDecayChannel.

In the former one we have the following check:

 if ( (G4MT_parent->GetParticleType() != "nucleus") &&
      (sumofdaughtermass > parentmass + 5*widthMass) ){
  // !!! illegal mass  !!!

the widthMass is the width of the parent particle. If we enter here the numbers for our problematic decay we get:

0.93956536 + 0.89581 < 1.8 + 5 * 300

so, the condition is satisfied and the decay is allowed. 
Btw, why the width of the daughters is not taken into account here??

OK, now we continue and we go to G4VDecayChannel where there is the check:

 daughtermomentum = Pmx(parentmass,daughtermass[0],daughtermass[1]);
 if (daughtermomentum <0.0) {
#ifdef G4VERBOSE
   if (GetVerboseLevel()>0) {
     G4cout << "G4PhaseSpaceDecayChannel::TwoBodyDecayIt " 
            << "sum of daughter mass is larger than parent mass" << G4endl;
     G4cout << "parent :" << G4MT_parent->GetParticleName() << "  " << current_parent_mass/GeV << G4endl;
     G4cout << "daughter 1 :" << G4MT_daughters[0]->GetParticleName() << "  " << daughtermass[0]/GeV << G4endl;
     G4cout << "daughter 2:" << G4MT_daughters[1]->GetParticleName() << "  " << daughtermass[1]/GeV << G4endl;
   }
#endif
   G4Exception("G4PhaseSpaceDecayChannel::TwoBodyDecayIt",
               "PART112", JustWarning,
               "Can not create decay products: sum of daughter mass is larger than parent mass”);



which does NOT take into account any width, only masses, and therefore it fails because 0.93956536 + 0.89581 > 1.8

So, this decay NEVER goes through, and the products of the decay are always empty.
The same happens for some other decays as well.

The simplest fix, is to make the second check consistent with the first one, i.e. to add the 5*widthMass, but actually all the widths should be taken into account.
Comment 1 Witold.Pokorski 2014-05-20 16:05:26 CEST
One typing mistake... The line:

"OK, now we continue and we go to G4VDecayChannel where there is the check:"

should read:

"OK, now we continue and we go to G4PhaseSpaceDecayChannel where there is the check:"
Comment 2 kurasige 2014-05-25 07:49:50 CEST
  Thanks to reporting problems.

  This is not a bug but I agree that we need to clarify how to treat mass and width in decay process.
 
  In G4VDecayChannel::FillDaughters(), sum of daughter mass is calculated and checked that it is consistent with parent mass within 5 times parent mass width. 

  In G4PhaseSpaceDecayChannel::DecayIt, momentum of daughter particles (in center of mass frame) are calculated according to given parent 'dynamical' mass (mass in DynamicParticle). In other words, the decay process can not modify the parent mass at all in order to keep energy momentum conservation.
 
 In cases such as anti_lambda(1800) decay, you need to give dynamical mass larger than PDG mass to the parent track when the parent track is created.

  As for treatment of daughter mass, I start to implement new G4PhaseSpaceDecayChannel, where daughter dynamic masses are assigned by considering parent mass and daughter mass width.
Comment 3 kurasige 2014-05-31 09:11:46 CEST
G4VDecayChannel and G4PhaseSpaceDecayChannel are modified to taking into account width of daughters particles. The dynamical mass of daughters particles will be given by G4PhaseSpaceDecayChannel so that sum of mass of daughter particles is less than the mass of the parent particle. The mass is calculated according to Breit-Wigner distribution (G4VDecayChannel::DynamicalMass method) within PDGMass +- PDGWidth.
(this mass range is defined by G4VDecalChannel::rangeMass)

Thanks to this modification, anti_lambda(1800) can decay into anti_neutron +  k_star0 even if the mass of anti_lambda is equal to PDGMass.  

The modification is included in the new tag of particle-V10-00-11