Changeset 1315 for trunk/source/processes/hadronic/models/chiral_inv_phase_space/processes/src/G4QLowEnergy.cc

Timestamp:

Jun 18, 2010, 11:42:07 AM (14 years ago)

Author:

garnier

Message:

update geant4-09-04-beta-cand-01 interfaces-V09-03-09 vis-V09-03-08

File:

• trunk/source/processes/hadronic/models/chiral_inv_phase_space/processes/src/G4QLowEnergy.cc (modified) (6 diffs)

trunk/source/processes/hadronic/models/chiral_inv_phase_space/processes/src/G4QLowEnergy.cc

@@ -24 +24 @@
 // ********************************************************************
 //
-// $Id: G4QLowEnergy.cc,v 1.1 2009/11/17 10:36:55 mkossov Exp $
-// GEANT4 tag $Name: geant4-09-03 $
+// $Id: G4QLowEnergy.cc,v 1.4 2010/06/07 15:19:55 mkossov Exp $
+// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
 //
 //      ---------------- G4QLowEnergy class -----------------
@@ -548 +548 @@
         if(!tD || !pD)                // Quasi-Elastic: B+A->(B-1)+N+A || ->B+N+(A-1)
         {
-          G4VQCrossSection* CSmanager=G4QElasticCrossSection::GetPointer();
+          G4VQCrossSection* PCSmanager=G4QProtonElasticCrossSection::GetPointer();
+          G4VQCrossSection* NCSmanager=G4QNeutronElasticCrossSection::GetPointer();
           G4int pPDG=2112;            // Proto of the nucleon PDG (proton)
           G4double prM =mNeut;        // Proto of the nucleon mass
@@ -587 +588 @@
               else break;             // Break the while loop
             }
-            G4double xSec=CSmanager->GetCrossSection(false, Mom, tZ, tN, pPDG);
+            G4double xSec=0.;
+            if(pPDG==2212) xSec=PCSmanager->GetCrossSection(false, Mom, tZ, tN, pPDG);
+            else           xSec=NCSmanager->GetCrossSection(false, Mom, tZ, tN, pPDG);
             if( xSec <= 0. ) break;   // Break the while loop
-            G4double mint=CSmanager->GetExchangeT(tZ,tN,pPDG); // functional randomized -t
-            G4double cost=1.-mint/CSmanager->GetHMaxT();       // cos(theta) in CMS
+            G4double mint=0.;                 // Prototype of functional randomized -t
+            if(pPDG==2212) mint=PCSmanager->GetExchangeT(tZ,tN,pPDG); // randomized -t
+            else           mint=NCSmanager->GetExchangeT(tZ,tN,pPDG); // randomized -t
+            G4double maxt=0.;                           // Prototype of maximum -t
+            if(pPDG==2212) maxt=PCSmanager->GetHMaxT(); // maximum -t
+            else           maxt=NCSmanager->GetHMaxT(); // maximum -t
+            G4double cost=1.-mint/maxt;       // cos(theta) in CMS
             if(cost>1. || cost<-1.) break; // Break the while loop
             G4LorentzVector reco4M=G4LorentzVector(0.,0.,0.,tgM); // 4mom of recoil target
@@ -729 +737 @@
               rNuc4M=G4LorentzVector(fm,rNE);
               G4LorentzVector qfN4M=targ4M - rNuc4M;// 4Mom of the quasi-free nucleon in LS
-              com4M += qfN4M;         // Calculate Total 4Mom for NA scattering
+              com4M += qfN4M;                 // Calculate Total 4Mom for NA scattering
               G4ThreeVector boostV=proj4M.vect()/proj4M.e();
               qfN4M.boost(-boostV);
-              G4double tNE = qfN4M.e(); // Energy of the QF nucleon in LS
+              G4double tNE = qfN4M.e();       // Energy of the QF nucleon in LS
               if(rNE >= prM) Mom = std::sqrt(tNE*tNE-prM2); // Mom(s) fake value
-              else break;             // Break the while loop
-            }
-            G4double xSec=CSmanager->GetCrossSection(false, Mom, pZ, pN, pPDG);
-            if( xSec <= 0. ) break;   // Break the while loop
-            G4double mint=CSmanager->GetExchangeT(pZ,pN,pPDG); // functional randomized -t
-            G4double cost=1.-mint/CSmanager->GetHMaxT();       // cos(theta) in CMS
-            if(cost>1. || cost<-1.) break; // Break the while loop
+              else break;                     // Break the while loop
+            }
+            G4double xSec=0.;
+            if(pPDG==2212) xSec=PCSmanager->GetCrossSection(false, Mom, tZ, tN, pPDG);
+            else           xSec=NCSmanager->GetCrossSection(false, Mom, tZ, tN, pPDG);
+            if( xSec <= 0. ) break;           // Break the while loop
+            G4double mint=0.;                 // Prototype of functional randomized -t
+            if(pPDG==2212) mint=PCSmanager->GetExchangeT(tZ,tN,pPDG); // randomized -t
+            else           mint=NCSmanager->GetExchangeT(tZ,tN,pPDG); // randomized -t
+            G4double maxt=0.;                           // Prototype of maximum -t
+            if(pPDG==2212) maxt=PCSmanager->GetHMaxT(); // maximum -t
+            else           maxt=NCSmanager->GetHMaxT(); // maximum -t
+            G4double cost=1.-mint/maxt;                 // cos(theta) in CMS
+            if(cost>1. || cost<-1.) break;    // Break the while loop
             G4LorentzVector reco4M=G4LorentzVector(0.,0.,0.,prM); // 4mom of recoil target
             G4LorentzVector scat4M=G4LorentzVector(0.,0.,0.,rNM); // 4mom of scattered N
@@ -747 +762 @@
             {
               G4cout<<"G4QLE::Tt="<<com4M.m()<<",p="<<prM<<",r="<<rNM<<",c="<<cost<<G4endl;
-              break;                  // Break the while loop
+              break;                          // Break the while loop
             }
             G4Track* targSpect = 0;
             G4Track* aNucTrack = 0;
-            if(tB > tD)               // Fill the residual nucleus
+            if(tB > tD)                       // Fill the residual nucleus
             {
               G4int rZ=tZ-tC;
@@ -928 +943 @@
             G4cout<<"G4QLE::PSDI: rpA=1, rpZ"<<rpZ<<G4endl;
 #endif
+          }
+          else if(rpA==2 && rpZ==0)            // nn decay
+          {
+            theDefinition = aNeutron;
+            G4LorentzVector f4M=G4LorentzVector(0.,0.,0.,mNeut); // 4mom of 1st neutron
+            G4LorentzVector s4M=G4LorentzVector(0.,0.,0.,mNeut); // 4mom of 2nd neutron
+#ifdef pdebug
+            G4cout<<"G4QLE::CPS->n+n,nn="<<rp4M.m()<<" >? 2*MNeutron="<<2*mNeutron<<G4endl;
+#endif
+            if(!G4QHadron(rp4M).DecayIn2(f4M, s4M))
+            {
+              G4cout<<"*W*G4QLE::CPS->n+n,t="<<rp4M.m()<<" >? 2*Neutron="<<2*mAlph<<G4endl;
+            }
+            G4DynamicParticle* pNeu = new G4DynamicParticle(theDefinition, f4M);
+            aFraPTrack = new G4Track(pNeu, localtime, position );
+            aFraPTrack->SetWeight(weight);                    //    weighted
+            aFraPTrack->SetTouchableHandle(trTouchable);
+            tt4M-=f4M;
+#ifdef edebug
+            totBaN-=2;
+            tch4M -=f4M;
+            G4cout<<">>G4QLEn::PSDI:n,tZ="<<totChg<<",tB="<<totBaN<<",t4M="<<tch4M<<G4endl;
+#endif
+#ifdef pdebug
+            G4cout<<"G4QLowEn::PSDI:-->ProjSpectA4M="<<f4M<<G4endl;
+#endif
+            ++nSec;
+            rp4M=s4M;
+          }
+          else if(rpA>2 && rpZ==0)            // Z=0 decay
+          {
+            theDefinition = aNeutron;
+            G4LorentzVector f4M=rp4M/rpA;     // 4mom of 1st neutron
+#ifdef pdebug
+            G4cout<<"G4QLE::CPS->Nn,M="<<rp4M.m()<<" >? N*MNeutron="<<rpA*mNeutron<<G4endl;
+#endif
+            for(G4int it=1; it<rpA; ++it)     // Fill (N-1) neutrons to output
+            {
+              G4DynamicParticle* pNeu = new G4DynamicParticle(theDefinition, f4M);
+              G4Track* aNTrack = new G4Track(pNeu, localtime, position );
+              aNTrack->SetWeight(weight);                    //    weighted
+              aNTrack->SetTouchableHandle(trTouchable);
+              result.push_back(aNTrack);
+            }
+            G4int nesc = rpA-1;
+            tt4M-=f4M*nesc;
+#ifdef edebug
+            totBaN-=nesc;
+            tch4M -=f4M*nesc;
+            G4cout<<">G4QLEn::PSDI:Nn,tZ="<<totChg<<",tB="<<totBaN<<",t4M="<<tch4M<<G4endl;
+#endif
+#ifdef pdebug
+            G4cout<<"G4QLowEn::PSDI:-->ProjSpectA4M="<<f4M<<G4endl;
+#endif
+            nSec+=nesc;
+            rp4M=f4M;
           }
           else if(rpA==8 && rpZ==4)            // Be8 decay