FAWP www1.gantep.tr/~hep/local/FAWP/ Fast Analysis With Pythia

1. FAWP http://www1.gantep.edu.tr/~hep/local/FAWP/ Fast Analysis With Pythia These notes will be updated as FAWP evolves Current version: 2010-11-24

2. What’s new? As of 06/11/2010 Add truth vertex vectors vx[], vy[], vz[] myNewParticle() now includes a list id (see later) The path is now /home/hep/fawp/ As of 10/11/2010: Vector mid[] is now renamed as mtrk[] It is the mother’s track number. id in Pythia is the pdg! Added vectors mmtrk[] and mmpdg[] The mother’s mother track number and pdg. As of 22/11/2010: A track’s properties is now a struct! This changes how you access a variable (see later). What you need to do: * Replace your copy of fawp.hand fawp.cpp/cc with the new version. * Replace your copy of data files with the new ones! * Fix your code for the path, myNewParticle()and mtrk[] changes.

3. Introduction FAWP provides a fast environment for looping over Pythia final state objects. With FAWP we can: - Investigate basic properties of LHC collisions - Study basic analysis techniques - Develop analysis techniques - Investigating event rates and perform trial analyses.

4. FAWP data files [Pythia] FAWP data files are generated using the Pythia event generator. FAWP algorithms run on data files, i.e. FAWP does not run Pythia, hence FAWP is very fast. Currently one file is available: minbias_10M_S3556767.dat (6.5 GB) 10 million minimum-bias events p-p collisions at 7 TeVcenter of mass. (other event types can be added where required) pythia.init( 2212, 2212, 7000.0); pythia.readString("SoftQCD:minBias = on"); pythia.readString("Random:setSeed = on"); pythia.readString("Random:seed = 3556767");

5. Data file - particle selection [Pythia] FAWP data files contain final state charged tracks and photons subject to the following cuts (designed to reduce file size): - pseudo-rapidity range || < 2.5 - charged transverse momentum > 1 GeV/c - photon energy > 1 GeV Too hard? if ( pythia.event[i].isFinal() && abs(pythia.event[i].eta()) < 2.5 && ( ( pythia.event[i].isCharged() && pythia.event[i].pT()>1.0 ) || ( pythia.event[i].id() == 22 && pythia.event[i].e() >1.0 ) ) ) { write the particle information }

6. Particle information [Pythia] The data file contains only generator-level information, one particle per line: pythia.event[i].px() Particle pythia.event[i].py() three- pythia.event[i].pz() momentum (GeV) pythia.event[i].id() Pdg code pythia.event[i].mother1() Mother track number pythia.event[pythia.event[i].mother1()].id() and pdg code and again for the mother of the mother. pythia.event[i].charge() Particle charge (e) pythia.event[i].xProd() Particle pythia.event[i].yProd() production pythia.event[i].zProd() vertex (mm)

7. Particle information (reconstructed level) After FAWP reads generator-level information from the data file, reconstructed-level information is constructed by applying a very simple model of the ATLAS detector resolution effects (see smearTrack.h): For charged tracks, transverse momentum (GeV/c) is Gaussian smeared with pTgiven by: pT / pT = 0.05% pT 1% For photons, energy (GeV) is Gaussian smeared with Egiven by : E / E = 10% / E1/2 1% Vertex resolution is not implemented yet!

8. FAWP particle vectors Particle information is stored in a vector: vector<Particle> tracks; An event of particles can be looped over as follows: for (inti=0; i<nTracks; i++) { } The i’th particle of vector tracks has the following information: tracks[i].tpxTruth-level x-momenta (double) tracks[i].tpyTruth-level y-momenta (double) tracks[i].tpzTruth-level z-momenta (double) tracks[i].pxReconstructed-level x-momenta (double) tracks[i].pyReconstructed-level y-momenta (double) tracks[i].pzReconstructed-level z-momenta (double) tracks[i].tvxTruth-level x-vertex (double) tracks[i].tvyTruth-level y-vertex (double) tracks[i].tvzTruth-level z-vertex (double) tracks[i].chgCharge (int) tracks[i].pdgPdg code (int) tracks[i].mtrkMother track number (int) tracks[i].mpdgMother pdg code (int) tracks[i].mmtrkMother’s mother track number (int) tracks[i].mmpdgMother’s mother pdg code (int) tracks[i].flagUser flag [initially zero] (int)

9. FAWP program structure FAWP is maintained on gul3.bim.gantep.edu.tr You can install it elsewhere by recreating the include directory: /home/hep/fawp/include/ fawp.hfawpStd.h randFlat.hrandGauss.hsmearTrack.hstdB.h You can take the main program from /home/hep/fawp/common/ compile-root For compiling fawp with root fawp.cc Example fawp with root fawp.cpp Example fawp without root You also need some data (maybe a subset of this): /home/hep/fawp/data/ minbias_10M_S3556767.dat

10. // fawp.cpp #include "/home/hep/fawp/include/fawp.h" using namespace fawp; int main() { openData(1); while (true) { getEvent(1000); if (eof) break; // EOF for (inti=0; i<nTracks; i++) { cout << tracks[i].tpx << endl; } // track loop } // event loop closeData(); } Basic structure of a FAWP program The data file is opened for reading. Event vectors are filled. (we request 1000 events here) Tracks are looped over, track vectors are processed. Not necessary, but you get some statistics

11. Other FAWP functions void myNewParticle( intmyid, double mypx, double mypy, double mypz, intmyflag=0, double mytpx=0.0, double mytpy=0.0, double mytpz=0.0 ) Use the myNewParticle()function to create lists of new derived particles. - myid (0 to 9) is used to identify the group that your derived particle belongs to (e.g. a list of 0s or D0s), this allows you to separately loop over different derived particle lists, and to create new particle lists while looping over a derived particle list. - A new (derived) particle has no pdg value or mother, charge and vertex. - The particle reconstructed 3-momenta should be given as parameters. - The user flag, and truth 3-momenta are optional and default to zero. - All FAWP functions work on derived particles created with myNewParticle() including derived particles derived from derived particles... The simplest form is myNewParticle(myid, mypx, mypy, mypz); To loop over your new particles in an event use the following loop structure: for (inti=myBegin[myid]; i<myEnd[myid]; i++) { } Note that if you wish to create a second derived particle list while looping over a first derived particle loop list then the second list myid must be greater than that of the first.

12. Q functions (work with track numbers) double Qp(inti) reconstructed momentum scaler for particle i double Qtp(inti) truth momentum scaler for particle i double QpT(inti) reconstructed transverse momentum scaler for particle i double QtpT(inti) truth transverse momentum scaler for particle i double Qeta(inti) reconstructed pseudo-rapidity for particle i double Qteta(inti) truth pseudo-rapidity for particle i double QtL(inti) truth vertex transverse distance from (vx=0,vy=0) for particle i double QtdL(inti, inj j) truth vertex transverse distance between particle i and j double Qm2(int i1, int i2, double m1=m_pi, double m2=m_pi) reconstructed invariant mass of particles i1 and i2. double Qtm2(int i1, int i2, double m1=m_pi, double m2=m_pi) truth invariant mass of particles i1 and i2. The third and fourth parameters are the hypotheses for the masses of particles i1 and i2 respectively; the default is m_pi. double Qm3(int i1, int i2, int i3, double m1=m_pi, double m2=m_pi, double m3=m_pi) double Qtm3(int i1, int i2, int i3, double m1=m_pi, double m2=m_pi, double m3=m_pi) as for Qm2 and Qtm2 but for three particles.

13. namespace fawpstd namespace fawpstd { inline double pow2(const double& x) {return x*x;} // particle properties directly from Pythia. // Masses are in GeV, mean lifelength in mm. const intid_e = 11; const double m_e = 0.0005110; const intid_mu = 13; const double m_mu = 0.10566, ct_mu = 658.654e+03; const intid_pi = 211; const double m_pi = 0.13957, ct_pi = 7.80450e+03; and so on for k, p, bs0, jpsi, phi, .... (you can request more to be added) }

14. Example - Plot the invariant mass of the decay X -> +-o // form pi0 candidates for ( i=0; i<nTracks-1; i++ ) { if ( track[i].chg==0 ) { for ( j=i+1; j<nTracks; j++ ) { if ( track[j].chg==0 ) { double qm2 = Qm2(i,j,0.0,0.0); h1->Fill( qm2 ); if ( qm2>0.110 && qm2<0.160 ) myNewParticle(0,pxi+pxj, pyi+pyj, pzi+pzj); }} }} h1 pyi=tracks[i].py for ( i=0; i<nTracks; i++ ) { if ( track[i].chg>0 ) { for ( j=0; j<nTracks; j++ ) { if ( track[j].chg<0 ) { for ( k=myBegin[0]; k<myEnd[0]; k++ ) h2->Fill( Qm3(i,j,k,m_pi,m_pi,0.135) ); }} }} h2 END