Subversion Repositories f9daq

Rev

Rev 266 | Go to most recent revision | Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
193 f9daq 1
#include <TClonesArray.h>
2
#include <TH1F.h>
3
#include <TStopwatch.h>
4
#include <TDatime.h>
5
#include <TString.h>
6
#include <TFile.h>
7
#include <TTree.h>
8
#include <TBranch.h>
9
#include <TClonesArray.h>
10
 
11
#include "BParticle.h"
12
#include "BEvent.h"
13
 
14
class Hdr{
15
public:
16
  int id;
17
  int len;
18
  int progress;
19
};
20
 
21
 
22
class Blab2 {
23
public:
24
UInt_t fNeve;
25
int fData;
26
TH1F *fH[100];
27
TClonesArray *fList[100];
28
 
29
Blab2();
30
~Blab2();
31
void Init();
32
void event();
33
void Process();
34
void h1d(const char *varname, const char *name, int nbins, double min, double max, int id );
35
int  selector(int pin, int charge, SIMPLEPID particlename,  int hid, int pout );
36
int  combiner(int id0 ,int id1 , SIMPLEPID particlename, double min, double max, int hid, int id );
37
int  fix_mass(int id);
38
void plist(int i);
39
 
40
 
41
ClassDef ( Blab2, 1 )
42
};
43
 
44
ClassImp(Blab2)
45
 
46
Blab2::Blab2():fNeve(0), fData(0){
47
 
48
 Process();
49
};
50
 
51
void Blab2::h1d(const char *varname, const char *name, int nbins, double min, double max, int id ){
52
   fH[id]= new TH1F(TString::Format("h%d",id), name, nbins, min, max);
53
}
54
 
55
int Blab2::combiner(int _p0, int _p1,SIMPLEPID pid, double min, double max, int hid, int _p2 ){
56
   // Loop over all the particles in both lists.
57
 if (_p0 < 0 ) _p0 =0;
58
 if (_p1 < 0 ) _p1 =0;
59
 
60
 
61
 fList[_p2]->Clear();
62
 int nprt=0;
63
 
64
 for(TIter next1(fList[_p0]);BParticle * p1 =(BParticle *) next1();) {
65
    // the second loop
66
 
67
   for(TIter next2 = (_p0!=_p1) ?  TIter(fList[_p1]): TIter(next1) ; BParticle * p2 =(BParticle *) next2();) {  
68
      if (p1==p2) continue;
69
      BParticle  p = *p1 + *p2; // Combine two particles into a new particle   
70
      if (p.InMassRange(min, max)){
71
            //SIMPLEPID pid = PION; // should be set to particlename
72
            p.SetPid(pid);
73
            TClonesArray& list = *fList[_p2];
74
            new (list[nprt++]) BParticle ( p ); // create a new entry in kslist list of particles
75
            if (hid>=0 && fH[hid]) fH[hid]->Fill(p.GetMass());
76
      }
77
 
78
   }
79
 
80
 }
81
 return _p2;
82
}
83
 
84
 
85
int Blab2::selector(int pin, int charge, SIMPLEPID type ,  int hid, int pout ){
86
 if (pin < 0 ) pin =0;
87
 
88
 
89
  fList[pout]->Clear();
90
  int nprt=0;
91
 
92
  for(TIter next(fList[pin]); BParticle * p =(BParticle *) next();) {
93
        if (p->charge()== charge && p->pid()== type ) {
94
          TClonesArray& list = *fList[pout];
95
          new (list[nprt++]) BParticle ( *p );
96
          if (hid>=0 && fH[hid]) fH[hid]->Fill(p->GetMass());
97
        }
98
  }      
99
   return pout;
100
}
101
 
102
 
103
int Blab2::fix_mass(int pin){
104
   if (pin < 0 ) pin =0;
105
   for(TIter next(fList[pin]); BParticle * p =(BParticle *) next();)  p->SetEnergyFromPid();
106
   return pin;
107
}
108
 
109
void Blab2::plist(int i){
110
  fList[i]= new TClonesArray( "BParticle", 500 );
111
}
112
Blab2::~Blab2(){};
113
 
114
 
115
void send_message(int id, TString msg, int progress){
116
static Hdr hdr;
117
 
118
   hdr.id = id;
119
   hdr.len= msg.Length();
120
   hdr.progress= progress;
121
   fwrite(&hdr,3*sizeof(int),1,stdout);
122
   fwrite(msg.Data(),hdr.len,1,stdout);
123
   fflush(stdout);
124
 
125
}
126
 
127
 
128
void Blab2::Process(){
129
 
130
for (int i=0;i<100;i++) fH[i]=NULL;
131
for (int i=0;i<100;i++) fList[i]=NULL;
132
 
133
Init();
134
 
135
TFile * f = new TFile(TString::Format("../../data/hadron-%d.root",fData)); // Open a data file
136
if(f->IsZombie()) {  send_message(0,TString::Format("File %d not found\n",fData), 0 );  return; }  
137
TTree * t =(TTree *) f-> Get( "T"); // Obtain a pointer to a tree of "event" data in the file
138
BEvent * mevent = new BEvent(); // Create a  "BEvent" object where data from the file will be loaded
139
TBranch * branch = t-> GetBranch( "BEvent"); // Obtain a branch for "BEvent" in the tree
140
branch-> SetAddress(&mevent); // Register created "BEvent" object to the branch
141
 
142
send_message(0, TString::Format("Number of Entries in the file %lld\n", t->GetEntries() ),0);
143
TStopwatch timer;
144
timer.Start();
145
int i=0;
146
int cNeve=TMath::Min(fNeve, (UInt_t) t-> GetEntries());
147
while (i<cNeve){
148
 t-> GetEntry(i); // Read the content of the event from the file
149
 fList[0]= mevent->GetParticleList();
150
 event();
151
 
152
 if (i%10000==0) send_message(2, TString::Format("Event %d\n",i), (100*i)/cNeve);
153
 mevent-> Clear();  // Clear the memory.
154
 for (int k=0;k<100;k++) if (fList[k]!=0) fList[k]->Clear();
155
 i++;
156
}
157
 
158
send_message(2, TString::Format("Number of events processed %d\n", i ),100);
159
for (int i=0;i<100;i++) if (fH[i]!=0) send_message(1,TBufferJSON::ConvertToJSON(fH[i]),100 );
160
 
161
TDatime d;
162
timer.Stop();
163
send_message(3, TString::Format("'%s', %d, %f, %f", d.AsSQLString(),i, timer.RealTime(), timer.CpuTime() ),100);
164
 
165
}
166
 
167