Realisation of AliGenReader that reads the old cwn event format.
[u/mrichter/AliRoot.git] / EVGEN / AliGenReaderCwn.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16
17 /*
18 $Log$
19 */
20
21 // Read the old ALICE event format based on CW-ntuples
22 // http://consult.cern.ch/alice/Internal_Notes/1995/32/abstract
23 // .cwn file have to be converted to .root using h2root
24 // Use SetFileName(file) to read from "file" 
25 // Author: andreas.morsch@cern.ch
26
27 #include <TFile.h>
28 #include <TTree.h>
29 #include <TParticle.h>
30
31 #include "AliGenReaderCwn.h"
32 #include "AliMC.h"
33 ClassImp(AliGenReaderCwn);
34
35
36 AliGenReaderCwn::AliGenReaderCwn() 
37 {
38 // Default constructor
39     fNcurrent   = 0;
40     fTreeNtuple = 0;
41 }
42
43 void AliGenReaderCwn::Init() 
44 {
45 //
46 // reset the existing file environment and open a new root file if
47 // the pointer to the Fluka tree is null
48     
49     TFile *pFile=0;
50     if (!pFile) {
51         pFile = new TFile(fFileName);
52         pFile->cd();
53         printf("\n I have opened %s file \n", fFileName);
54     }
55 // get the tree address in the Fluka boundary source file
56     fTreeNtuple = (TTree*)gDirectory->Get("h888");
57
58     TTree *h2=fTreeNtuple;
59 //Set branch addresses
60     h2->SetBranchAddress("Nihead",&fNihead);
61     h2->SetBranchAddress("Ihead",fIhead);
62     h2->SetBranchAddress("Nrhead",&fNrhead);
63     h2->SetBranchAddress("Rhead",fRhead);
64     h2->SetBranchAddress("Idpart",&fIdpart);
65     h2->SetBranchAddress("Theta",&fTheta);
66     h2->SetBranchAddress("Phi",&fPhi);
67     h2->SetBranchAddress("P",&fP);
68     h2->SetBranchAddress("E",&fE);
69 }
70
71 Int_t AliGenReaderCwn::NextEvent() 
72 {
73 // Read the next event  
74     Int_t nTracks;
75     fNparticle = 0;
76     TFile* pFile = fTreeNtuple->GetCurrentFile();
77     pFile->cd();
78
79     Int_t nentries = (Int_t) fTreeNtuple->GetEntries();
80     if (fNcurrent < nentries) {
81         Int_t nb = (Int_t)fTreeNtuple->GetEvent(fNcurrent);
82         fNcurrent++;
83         
84         Int_t i5=fIhead[4];
85         Int_t i6=fIhead[5];
86         if (i5==0) {
87             printf("\n This should never happen !\n");
88             nTracks = 0;
89         } else {
90             printf("\n Next event contains %d tracks! \n", i6);
91             nTracks = i6;
92         }    
93         fNparticleMax = nTracks;
94         return nTracks;
95     } else {
96         return 0;
97     }
98     return 0;
99 }
100
101 TParticle* AliGenReaderCwn::NextParticle() 
102 {
103 //
104 //  
105     Float_t prwn;
106     Float_t p[4];
107 // Read the next particle
108     if (fCode == kGEANT3) fIdpart=gMC->PDGFromId(fIdpart);
109     Double_t amass = TDatabasePDG::Instance()->GetParticle(fIdpart)->Mass();
110     if(fE<=amass) {
111         Warning("Generate","Particle %d  E = %f mass = %f %f %f \n",
112                 fIdpart,fE,amass, fPhi, fTheta);
113         prwn=0;
114     } else {
115         prwn=sqrt((fE+amass)*(fE-amass));
116     }
117
118     fTheta *= TMath::Pi()/180.;
119     fPhi    = (fPhi-180)*TMath::Pi()/180.;      
120     p[0] = prwn*TMath::Sin(fTheta)*TMath::Cos(fPhi);
121     p[1] = prwn*TMath::Sin(fTheta)*TMath::Sin(fPhi);      
122     p[2] = prwn*TMath::Cos(fTheta);
123     p[3] = fE;
124     TParticle* particle = new TParticle(fIdpart, 0, -1, -1, -1, -1, p[0], p[1], p[2], p[3], 
125                                         0., 0., 0., 0.);
126     Int_t nb = (Int_t)fTreeNtuple->GetEvent(fNcurrent);
127     fNcurrent++;
128     fNparticle++;
129     return particle;
130 }
131
132
133
134
135
136