WebSVN – f9daq – Diff – /belle2/masterclass/src/Blab2.cc

 #include <TDatime.h>
 #include <TString.h>
 #include <TFile.h>
 #include <TTree.h>
 #include <TBranch.h>
-#include <TClonesArray.h>
+#include <TBufferJSON.h>
 #include "BParticle.h"
 #include "BEvent.h"
 class Hdr{
 };
 class Blab2 {
 public:
+const char *names[11]={"photon", "electron", "pion", "muon", "kaon", "proton", "J/Psi", "D", "D*", "B", "Phi"};
 UInt_t fNeve;
 UInt_t fNfirst;
 UInt_t fPrint;
 int fData;
 TH1F *fH[100];
 };
 void Blab2::h1d(const char *varname, const char *name, int nbins, double min, double max, int id ){
    TString svar(varname);
-   TString axis[3]={"mass (GeV/c2)","momentum (GeV/c)","energy (GeV)"};
+   TString axis[]={"mass (GeV/c2)","momentum (GeV/c)","energy (GeV)","charge","identity","momentum (GeV/c)","momentum (GeV/c)","momentum (GeV/c)","momentum (GeV/c)"};
    fHtype[id] = 0;
    if (svar.Contains("GetMass"    )) fHtype[id]=0;
    if (svar.Contains("GetMomentum")) fHtype[id]=1;
    if (svar.Contains("GetEnergy"  )) fHtype[id]=2;
+   if (svar.Contains("GetCharge"  )) fHtype[id]=3;
+   if (svar.Contains("GetPid"     )) fHtype[id]=4;
+   if (svar.Contains("GetXMomentum")) fHtype[id]=5;
+   if (svar.Contains("GetYMomentum")) fHtype[id]=6;
+   if (svar.Contains("GetZMomentum")) fHtype[id]=7;
+   if (svar.Contains("GetTranverseMomentum")) fHtype[id]=8;
    fH[id]= new TH1F(TString::Format("h%d",id), TString::Format("%s;%s;N",name,axis[fHtype[id]].Data()), nbins, min, max);
+   if (fHtype[id]==4) {
+     fH[id]= new TH1F(TString::Format("h%d",id), TString::Format("%s;%s;N",name,axis[fHtype[id]].Data()), 11, 0, 11);
+     for (int i=0;i<11;i++) fH[id]->GetXaxis()->SetBinLabel(i+1,names[i]);
+   } else {
+     fH[id]= new TH1F(TString::Format("h%d",id), TString::Format("%s;%s;N",name,axis[fHtype[id]].Data()), nbins, min, max);
+   }
+}
 int Blab2::Fill(int hid, BParticle *p){
 if (hid>=0 && fH[hid]) {
       float val;
       switch (fHtype[hid]){
       case 0: val = p->GetMass(); break;
       case 1: val = p->GetMomentum(); break;
       case 2: val = p->e(); break;
+      case 3: val = p->charge(); break;
+      case 4: val = p->pid(); break;
+      case 5: val = p->px(); break;
+      case 6: val = p->py(); break;
+      case 7: val = p->pz(); break;
+      case 8: val = sqrt(p->px()*p->px()+p->py()*p->py()); break;
       default : val = 0 ; break;
+      }
    fH[hid]->Fill(val);
+}
  if (_p0 < 0 ) _p0 =0;
  if (_p1 < 0 ) _p1 =0;
  fList[_p2]->Clear();
  int nprt=0;
  for(TIter next1(fList[_p0]);BParticle * p1 =(BParticle *) next1();) {
     // the second loop
    // in the case the second parti
    for(TIter next2 = (_p0!=_p1 && same==0) ?  TIter(fList[_p1]): TIter(next1) ; BParticle * p2 =(BParticle *) next2();) {
       if (p1==p2) continue;     // do not use the same particle in the combinations
       BParticle  p = *p1 + *p2; // Combine two particles into a new particle
       if (p.InMassRange(min, max)){
+            Fill(hid, &p);
             p.SetPid(pid); // set PID to particlename to fix the particle mass
+            p.SetEnergyFromPid();
             TClonesArray& list = *fList[_p2];
             new (list[nprt++]) BParticle ( p ); // create a new entry in kslist list of particles
-            Fill(hid, &p);
+      }
+   }
+ }
+}
 int Blab2::selector(int pin, int charge, SIMPLEPID type ,  int hid, int pout ){
  if (pin < 0 ) pin =0;
   fList[pout]->Clear();
   int nprt=0;
   for(TIter next(fList[pin]); BParticle * p =(BParticle *) next();) {
+        if (p->charge()== charge || charge > 1){
-        if (p->charge()== charge && ( p->pid()== type || type == ALL )) {
+          if ( p->pid()== type || type == ALL ) {
-          TClonesArray& list = *fList[pout];
+            TClonesArray& list = *fList[pout];
-          new (list[nprt++]) BParticle ( *p );
+            new (list[nprt++]) BParticle ( *p );
-          Fill(hid, p);
+            Fill(hid, p);
+          }
+        }
+  }
    return pout;
 }
 int Blab2::fix_mass(int pin){
    if (pin < 0 ) pin =0;
    for(TIter next(fList[pin]); BParticle * p =(BParticle *) next();)  p->SetEnergyFromPid();
    return pin;
 }
 void Blab2::plist(int i){
   fList[i]= new TClonesArray( "BParticle", 500 );
+}
 Blab2::~Blab2(){};
 void send_message(int id, TString msg, int progress){
 static Hdr hdr;
    hdr.id = id;
    hdr.len= msg.Length();
    hdr.progress= progress;
    fwrite(&hdr,3*sizeof(int),1,stdout);
    fwrite(msg.Data(),hdr.len,1,stdout);
 for (int i=0;i<100;i++) fH[i]=NULL;
 for (int i=0;i<100;i++) fHtype[i]=0;
 for (int i=0;i<100;i++) fList[i]=NULL;
 Init();
 TFile * f = new TFile(TString::Format("../../data/hadron-%d.root",fData)); // Open a data file
 if(f->IsZombie()) {  send_message(0,TString::Format("File %d not found\n",fData), 0 );  return; }
 TTree * t =(TTree *) f-> Get( "T"); // Obtain a pointer to a tree of "event" data in the file
 BEvent * mevent = new BEvent(); // Create a  "BEvent" object where data from the file will be loaded
 TBranch * branch = t-> GetBranch( "BEvent"); // Obtain a branch for "BEvent" in the tree
 branch-> SetAddress(&mevent); // Register created "BEvent" object to the branch
 TH1F *fHnprt= new TH1F("h100", "Number of particles in the event;N particles;N events", 50, -0.5, 49.5);
-TH1F *fHid= new TH1F("h101", "Particle types;ID;N particles", 6, -0, 6);
-const char *names[12]={"PHOTON", "ELECTRON", "PION", "MUON", "KAON", "PROTON", "JPSI", "D", "DSTAR", "B","KS", "ALL" };
-for (int i=1;i<=6;i++) fHid->GetXaxis()->SetBinLabel(i,names[i-1]);
 send_message(0, TString::Format("<br>Number of Events in the file %lld<br>\n", t->GetEntries() ),0);
 TStopwatch timer;
 timer.Start();
 int nev  = 0;
  int nprt=0;
  if (nev>100) fPrint = 0; // disable particle prints for huge numer of events
  if (fPrint) sprintf(sList,"Primary particle list for Event %d<br/><table class='plist' ><tr><th>N<th>px(GeV/c)<th>py(GeV/c)<th>pz(GeV/c)<th>p(GeV/c)<th>Energy(GeV)<th>Charge<th>ID<th></tr>", i);
  for(TIter next(fList[0]); BParticle * p =(BParticle *) next();) {
    nprt++;
-   fHid->Fill(p->pid());
    if (fPrint) sprintf(sList,"%s<tr><td>%d<td>%g<td>%g<td>%g<td>%g<td>%g<td>%1.0f<td>%s</tr>",sList,nprt, p->px(),p->py(),p->pz(),p->GetMomentum(),p->e(), p->charge(),names[p->pid()] );
+ }
  fHnprt->Fill(nprt);
  if (fPrint) {
    sprintf(sList,"%s</table>",sList);
+}
 send_message(2, TString::Format("Number of events processed %d\n", i ),100);
 send_message(1,TBufferJSON::ConvertToJSON(fHnprt),100 );
-send_message(1,TBufferJSON::ConvertToJSON(fHid),100 );
 for (int i=0;i<100;i++) if (fH[i]!=0) send_message(1,TBufferJSON::ConvertToJSON(fH[i]),100 );
 TDatime d;
 timer.Stop();

Subversion Repositories f9daq

(root)/belle2/masterclass/src/Blab2.cc – Rev 268 → 272