source: trunk/source/processes/hadronic/models/binary_cascade/test/G4XToTest.cc@ 1228

Last change on this file since 1228 was 1199, checked in by garnier, 16 years ago

nvx fichiers dans CVS
File size: 8.8 KB
Line 
1//
2// ********************************************************************
3// * DISCLAIMER *
4// * *
5// * The following disclaimer summarizes all the specific disclaimers *
6// * of contributors to this software. The specific disclaimers,which *
7// * govern, are listed with their locations in: *
8// * http://cern.ch/geant4/license *
9// * *
10// * Neither the authors of this software system, nor their employing *
11// * institutes,nor the agencies providing financial support for this *
12// * work make any representation or warranty, express or implied, *
13// * regarding this software system or assume any liability for its *
14// * use. *
15// * *
16// * This code implementation is the intellectual property of the *
17// * GEANT4 collaboration. *
18// * By copying, distributing or modifying the Program (or any work *
19// * based on the Program) you indicate your acceptance of this *
20// * statement, and all its terms. *
21// ********************************************************************
22//
23//
24// -------------------------------------------------------------------
25// GEANT 4 class file --- Copyright CERN 1998
26// CERN Geneva Switzerland
27//
28//
29// File name: G4XResonanceTest.cc
30//
31// Author: Maria Grazia Pia (pia@genova.infn.it),
32//
33// Creation date: 15 June 1999
34//
35// Modifications:
36//
37// -------------------------------------------------------------------
38
39#include "globals.hh"
40
41#include "G4ios.hh"
42#include <fstream>
43#include <iomanip>
44#include <iostream>
45#include <assert.h>
46
47#include "CLHEP/Hist/TupleManager.h"
48#include "CLHEP/Hist/HBookFile.h"
49#include "CLHEP/Hist/Histogram.h"
50#include "CLHEP/Hist/Tuple.h"
51
52#include "Randomize.hh"
53
54#include "G4ThreeVector.hh"
55#include "G4LorentzVector.hh"
56#include "G4LorentzRotation.hh"
57
58#include "G4AntiProton.hh"
59#include "G4AntiNeutron.hh"
60#include "G4Proton.hh"
61#include "G4Neutron.hh"
62#include "G4PionPlus.hh"
63#include "G4PionMinus.hh"
64#include "G4PionZero.hh"
65#include "G4Gamma.hh"
66#include "G4MuonMinus.hh"
67#include "G4MuonPlus.hh"
68#include "G4KaonMinus.hh"
69#include "G4KaonPlus.hh"
70#include "G4NeutrinoMu.hh"
71#include "G4AntiNeutrinoMu.hh"
72#include "G4Lambda.hh"
73
74#include "G4ParticleDefinition.hh"
75
76#include "G4VDecayChannel.hh"
77#include "G4DecayTable.hh"
78
79#include "G4KineticTrack.hh"
80#include "G4KineticTrackVector.hh"
81
82#include "G4ParticleTable.hh"
83#include "G4VShortLivedParticle.hh"
84#include "G4ShortLivedConstructor.hh"
85#include "G4ShortLivedTable.hh"
86
87#include "G4VCrossSectionSource.hh"
88#include "G4XAqmElastic.hh"
89#include "G4XNNTotalLowE.hh"
90#include "G4XPDGTotal.hh"
91#include "G4XNNTotal.hh"
92#include "G4XPDGElastic.hh"
93
94#include "G4VXResonanceTable.hh"
95
96#include "G4XResonance.hh"
97#include "G4XDeltaDeltaToNN.hh"
98#include "G4XDeltaDeltaToNDelta.hh"
99#include "G4XNDeltaToNN.hh"
100#include "G4XNDeltaToDeltaDelta.hh"
101
102int main()
103{
104 // MGP ---- HBOOK initialization
105 HepTupleManager* hbookManager;
106 hbookManager = new HBookFile("xresonance.hbook", 58);
107 assert (hbookManager != 0);
108
109 // MGP ---- Book histograms
110
111 HepHistogram* hSigma;
112 hSigma = hbookManager->histogram("CrossSection", 100,0.,10.);
113 assert (hSigma != 0);
114
115 // MGP ---- Book a ntuple
116 HepTuple* ntuple;
117 ntuple = hbookManager->ntuple("CrossSection ntuple");
118 assert (ntuple != 0);
119
120 G4ParticleDefinition* gamma = G4Gamma::GammaDefinition();
121
122 G4ParticleDefinition* proton = G4Proton::ProtonDefinition();
123 G4ParticleDefinition* antiProton = G4AntiProton::AntiProtonDefinition();
124 G4ParticleDefinition* neutron = G4Neutron::NeutronDefinition();
125 G4ParticleDefinition* antiNeutron = G4AntiNeutron::AntiNeutronDefinition();
126
127 G4ParticleDefinition* pionPlus = G4PionPlus::PionPlusDefinition();
128 G4ParticleDefinition* pionMinus = G4PionMinus::PionMinusDefinition();
129 G4ParticleDefinition* pionZero = G4PionZero::PionZeroDefinition();
130
131 G4ParticleDefinition* kaonPlus = G4KaonPlus::KaonPlusDefinition();
132 G4ParticleDefinition* kaonMinus = G4KaonMinus::KaonMinusDefinition();
133
134 G4ParticleDefinition* lambda = G4Lambda::LambdaDefinition();
135
136 G4ParticleDefinition* theMuonPlus = G4MuonPlus::MuonPlusDefinition();
137 G4ParticleDefinition* theMuonMinus = G4MuonMinus::MuonMinusDefinition();
138
139 G4ParticleDefinition* theNeutrinoMu = G4NeutrinoMu::NeutrinoMuDefinition();
140 G4ParticleDefinition* theAntiNeutrinoMu = G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
141
142 G4ShortLivedConstructor ShortLived;
143 ShortLived.ConstructParticle();
144
145
146 // N1440
147 G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
148 G4ParticleDefinition* N1440Zero = particleTable->FindParticle("N(1440)0");
149 G4cout << N1440Zero->GetParticleName() << " created, type is "
150 << N1440Zero->GetParticleType() << G4endl;
151 G4ParticleDefinition* N1440Plus = particleTable->FindParticle("N(1440)+");
152 G4cout << N1440Plus->GetParticleName() << " created, type is "
153 << N1440Plus->GetParticleType() << G4endl;
154
155 // Delta
156 G4ParticleDefinition* deltaPlus = particleTable->FindParticle("delta+");
157 G4cout << deltaPlus->GetParticleName() << " created, type is "
158 << deltaPlus->GetParticleType() << G4endl;
159 G4String name = N1440Plus->GetParticleName();
160
161 // Select track 1
162
163 G4ParticleDefinition* definition1 = proton;
164 G4KineticTrack trk1;
165 G4double mass1 = definition1->GetPDGMass();
166
167 // Formation time
168 G4double time1 = 0.;
169 // G4cout << "Formation time (ns)" << G4endl;
170 // G4cin >> time1;
171 // if (time1 < 0.) time1 = 0.;
172 // time1 = time1 * ns;
173
174 // Initial position
175 G4double x1 = 0.;
176 G4double y1 = 0.;
177 G4double z1 = 0.;
178 // G4cout <<"Initial position (cm) " << G4endl;
179 // G4cin >> x1 >> y1 >> z1;
180 G4ThreeVector position1(x1*cm, y1*cm, z1*cm);
181
182 // Momentum
183 G4double pzMin1;
184
185 G4double pzMax1;
186 G4cout << "Pz min and max (GeV)" << G4endl;
187 G4cin >> pzMin1 >> pzMax1;
188 pzMin1 = pzMin1 * GeV;
189 pzMax1 = pzMax1 * GeV;
190
191 // Select track 2
192
193 G4ParticleDefinition* definition2 = deltaPlus;
194 G4KineticTrack trk2;
195
196 G4double mass2 = definition2->GetPDGMass();
197
198 // Formation time
199 G4double time2 = 0.;
200 // G4cout << "Formation time (ns)" << G4endl;
201 // G4cin >> time2;
202 // if (time2 < 0.) time2 = 0.;
203 // time2 = time2 * ns;
204
205 // Initial position
206 G4double x2 = 0.;
207 G4double y2 = 0.;
208 G4double z2 = 0.;
209 // G4cout <<"Initial position (cm) " << G4endl;
210 // G4cin >> x2 >> y2 >> z2;
211 G4ThreeVector position2(x2*cm, y2*cm, z2*cm);
212
213 // Momentum
214 G4double pzMin2;
215 G4double pzMax2;
216 G4cout << "Pz min and max (GeV)" << G4endl;
217 G4cin >> pzMin2 >> pzMax2;
218 pzMin2 = pzMin2 * GeV;
219 pzMax2 = pzMax2 * GeV;
220
221
222 // Create the ResonanceCrossSection
223
224 // G4VCrossSectionSource* source = new G4XDeltaDeltaToNDelta(definition1,definition2,
225 // deltaPlus,N1440Plus);
226
227 // G4VCrossSectionSource* source = new G4XDeltaDeltaToNN(definition1,definition1,
228 // proton,proton);
229
230 G4VCrossSectionSource* source = new G4XNDeltaToNN(proton,definition2,
231 proton,proton);
232
233 G4cout << "***** " << source->Name() << "***** created " << G4endl;
234 source->Print();
235
236 // Number of iterations
237 G4int n = 0;
238 G4cout << "---- Number of iterations ----" << G4endl;
239 G4cin >> n;
240 if (n < 1) n = 1;
241
242 G4double step1 = (pzMax1 - pzMin1) / n;
243 G4double step2 = (pzMax2 - pzMin2) / n;
244
245 G4cout << "pz min max step "
246 << pzMin1 / GeV << ", "
247 << pzMax1 / GeV << ", "
248 << step1 / GeV << ", " << G4endl
249 << "pz min max step "
250 << pzMin2 / GeV << ", "
251 << pzMax2 / GeV << ", "
252 << step2 / GeV << G4endl;
253
254 G4int i;
255 for (i=0; i<n; i++)
256 {
257 G4double pStep1 = pzMin1 + step1 * i;
258 G4ThreeVector p1(0.,0.,pStep1);
259 G4double e1 = sqrt(mass1*mass1 + p1.mag()*p1.mag());
260 G4LorentzVector p41(p1,e1);
261
262 G4KineticTrack trk1(definition1,time1,position1,p41);
263
264 G4double pStep2 = -(pzMin2 + step2 * i);
265 G4ThreeVector p2(0.,0.,pStep2);
266 G4double e2 = sqrt(mass2*mass2 + p2.mag()*p2.mag());
267 G4LorentzVector p42(p2,e2);
268 G4KineticTrack trk2(definition2,time2,position2,p42);
269
270 G4double sqrtS = (trk1.Get4Momentum() + trk2.Get4Momentum()).mag();
271
272 G4double sigma = source->CrossSection(trk1,trk2);
273
274 G4bool dummy = false;
275
276 G4cout << definition1->GetParticleName() << " " << definition2->GetParticleName()
277 << " cross section = " << sigma / millibarn << G4endl;
278
279 // ntuple->column("ecm",sqrtS/GeV);
280 // ntuple->column("sigma",sigma/millibarn);
281 // ntuple->dumpData();
282 }
283
284
285 hbookManager->write();
286
287 delete source;
288
289 return EXIT_SUCCESS;
290}
291
292
Note: See TracBrowser for help on using the repository browser.