Subversion Repositories f9daq

Rev

Rev 266 | Go to most recent revision | Blame | Compare with Previous | Last modification | View Log | RSS feed

  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.  
  168.