Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/inclxx/incl_physics/src/G4INCLPiNToMultiPionsChannel.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  1 //
  2 // ********************************************************************
  3 // * License and Disclaimer                                           *
  4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.                             *
 10 // *                                                                  *
 11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                                                  *
 18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // ********************************************************************
 25 //
 26 // INCL++ intra-nuclear cascade model
 27 // Alain Boudard, CEA-Saclay, France
 28 // Joseph Cugnon, University of Liege, Belgium
 29 // Jean-Christophe David, CEA-Saclay, France
 30 // Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland
 31 // Sylvie Leray, CEA-Saclay, France
 32 // Davide Mancusi, CEA-Saclay, France
 33 //
 34 #define INCLXX_IN_GEANT4_MODE 1
 35 
 36 #include "globals.hh"
 37 
 38 #include "G4INCLPiNToMultiPionsChannel.hh"
 39 #include "G4INCLKinematicsUtils.hh"
 40 #include "G4INCLBinaryCollisionAvatar.hh"
 41 #include "G4INCLRandom.hh"
 42 #include "G4INCLGlobals.hh"
 43 #include "G4INCLLogger.hh"
 44 #include <algorithm>
 45 #include "G4INCLPhaseSpaceGenerator.hh"
 46 
 47 namespace G4INCL {
 48 
 49   const G4double PiNToMultiPionsChannel::angularSlope = 15.;
 50 
 51   PiNToMultiPionsChannel::PiNToMultiPionsChannel(const G4int npi, Particle *p1, Particle *p2)
 52     : npion(npi),
 53     ind2(0),
 54     particle1(p1),
 55     particle2(p2)
 56   {
 57     std::fill(isosp, isosp+4, 0);
 58   }
 59 
 60   PiNToMultiPionsChannel::~PiNToMultiPionsChannel(){
 61 
 62   }
 63 
 64   void PiNToMultiPionsChannel::fillFinalState(FinalState *fs) {
 65 
 66 // assert(npion > 1 && npion < 5);
 67 
 68     Particle * nucleon;
 69     Particle * pion;
 70     if(particle1->isNucleon()) {
 71       nucleon = particle1;
 72       pion = particle2;
 73     } else {
 74       nucleon = particle2;
 75       pion = particle1;
 76     }
 77 #ifdef INCLXX_IN_GEANT4_MODE
 78     // Erase the parent resonance information of the nucleon and pion
 79     nucleon->setParentResonancePDGCode(0);
 80     nucleon->setParentResonanceID(0);
 81     pion->setParentResonancePDGCode(0);
 82     pion->setParentResonanceID(0);
 83 #endif
 84       G4int ipi=ParticleTable::getIsospin(pion->getType());
 85       ind2=ParticleTable::getIsospin(nucleon->getType());
 86 
 87       ParticleList list;
 88       list.push_back(nucleon);
 89       list.push_back(pion);
 90       fs->addModifiedParticle(nucleon);
 91       fs->addModifiedParticle(pion);
 92 
 93       isospinRepartition(ipi);
 94 
 95       const ParticleType tn=ParticleTable::getNucleonType(ind2);
 96       nucleon->setType(tn);
 97       ParticleType pionType=ParticleTable::getPionType(isosp[0]);
 98       pion->setType(pionType);
 99       const ThreeVector &rcolpion = pion->getPosition();
100       const ThreeVector zero;
101       for(G4int i=1; i<npion; ++i) {
102         pionType=ParticleTable::getPionType(isosp[i]);
103         Particle *newPion = new Particle(pionType,zero,rcolpion);
104         newPion->setType(pionType);
105         list.push_back(newPion);
106         fs->addCreatedParticle(newPion);
107       }
108 
109       const G4double sqrtS = KinematicsUtils::totalEnergyInCM(nucleon, pion);
110       PhaseSpaceGenerator::generateBiased(sqrtS, list, 0, angularSlope);
111 
112   }
113 
114     void PiNToMultiPionsChannel::isospinRepartition(G4int ipi) {
115         G4double rjcd=Random::shoot();
116         const G4int itot=ipi*ind2;
117 
118         isosp[1]=ipi;
119         if (npion != 3) {
120             if (npion == 4){
121               const G4double r2 = Random::shoot();
122                 if (r2*3. > 2.) {
123                     isosp[2]= 0;
124                     isosp[3]= 0;
125                 }
126                 else {
127                     isosp[2]=  2;
128                     isosp[3]= -2;
129                 }
130             }
131 
132             if (itot == 2) {
133                 // CAS PI+ P ET PI- n
134                 rjcd *= 5.;
135                 if (rjcd > 3.) {
136                     //      PI+ PI+ N ET PI- PI- P
137                     isosp[0]=2*ind2;
138                     isosp[1]=ipi;
139                     ind2=-ind2;
140                 }
141                 else {
142                     //      PI+ PI0 P ET PI- PI0 N
143                     isosp[0]=0;
144                     isosp[1]=ipi;
145                 }
146             }
147             else if (itot == 0) {
148                 // CAS PI0 P ET PI0 N
149                 rjcd *= 90.;
150                 if (rjcd > 13.) {
151                     if (rjcd > 52.) {
152                         //     PI+ PI0 N ET PI- PI0 P
153                         isosp[0]=2*ind2;
154                         isosp[1]=0;
155                         ind2=-ind2;
156                     }
157                     else {
158                         //     PI+ PI- P ET PI+ PI- N
159                         isosp[1]=-2;
160                         isosp[0]=2;
161                     }
162                 }
163                 else {
164                     //     PI0 PI0 P ET PI0 PI0 N
165                     isosp[0]=0;
166                     isosp[1]=0;
167                 }
168             }
169             else if (itot == -2) {
170                 // CAS PI- P ET PI+ N
171                 rjcd *= 45.;
172                 if (rjcd > 17.) {
173                     if (rjcd > 24.) {
174                         //     PI+ PI- N ET PI+ PI- P
175                         isosp[0]=2*ind2;
176                         ind2=-ind2;
177                     }
178                     else {
179                         //     PI0 PI0 N ET PI0 PI0 P
180                         isosp[0]=0;
181                         isosp[1]=0;
182                         ind2=-ind2;
183                     }
184                 }
185                 else
186                     //     PI- PI0 P ET PI+ PI0 N
187                     isosp[0]=0;
188             }
189         } // if (npion != 3)
190         else {
191             // PI N -> PI PI PI N
192             //          IF (IPI*IND2) 20,21,22
193             if (itot == -2) {
194                 // CAS PI- P ET PI+ N
195                 rjcd *= 135.;
196                 if (rjcd <= 28.) {
197                     //     PI0 PI0 PI0 N ET PI0 PI0 PI0 P
198                     isosp[0]=0;
199                     isosp[1]=0;
200                     isosp[2]=0;
201                     ind2=-ind2;
202                 }
203                 else {
204                     if (rjcd <= 84.) {
205                         //     PI+ PI- PI0 N ET PI- PI+ PI0 P
206                         isosp[0]=2*ind2;
207                         isosp[2]=0;
208                         ind2=-ind2;
209                     }
210                     else {
211                         if (rjcd <= 118.) {
212                             //     PI- PI- PI+ P ET PI- PI+ PI+ N
213                             isosp[0]=ipi;
214                             isosp[2]=-ipi;
215                         }
216                         else {
217                             //     PI- PI0 PI0 P ET PI0 PI0 PI+ N
218                             isosp[0]=0;
219                             isosp[2]=0;
220                         }
221                     }
222                 }
223             }
224             else if (itot == 0) {
225                 // CAS PI0 P ET PI0 N
226                 rjcd *= 270.;
227                 if (rjcd <= 39.) {
228                     //     PI0 PI0 PI0 P ET PI0 PI0 PI0 N
229                     isosp[0]=0;
230                     isosp[2]=0;
231                 }
232                 else {
233                     if (rjcd <= 156.) {
234                         //     PI+ PI- PI0 P ET PI- PI+ PI0 N
235                         isosp[0]=2;
236                         isosp[2]=-2;
237                     }
238                     else {
239                         if (rjcd <= 194.) {
240                             //     PI+ PI0 PI0 N ET PI0 PI0 PI- P
241                             isosp[0]=0;
242                             isosp[2]=2*ind2;
243                             ind2=-ind2;
244                         }
245                         else {
246                             //     PI- PI+ PI+ N ET PI- PI- PI+ P
247                             isosp[0]=2*ind2;
248                             isosp[1]=2*ind2;
249                             isosp[2]=-2*ind2;
250                             ind2=-ind2;
251                         }
252                     }
253                 }
254             }
255             else if (itot == 2) {
256                 // CAS PI+ P ET PI- N
257                 rjcd *= 5.;
258                 if (rjcd <= 2.) {
259                     //     PI+ P PI0 PI0 ET PI- N PI0 PI0
260                     isosp[0]=0;
261                     isosp[2]=0;
262                 }
263                 else {
264                     if (rjcd <= 3.) {
265                         //     PI+ P PI+ PI- ET PI- N PI+ PI-
266                         isosp[0]=-2;
267                         isosp[2]=2;
268                     }
269                     else {
270                         //     PI+ N PI+ PI0 ET PI- P PI0 PI-
271                         isosp[0]=2*ind2;
272                         isosp[2]=0;
273                         ind2=-ind2;
274                     }
275                 }
276             }
277         }
278 
279         std::shuffle(isosp,isosp+npion,Random::getAdapter()); // isospin randomly distributed
280     }
281 
282 }
283