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

 #include <TString.h>
 #include <TFile.h>
 #include <TTree.h>
 #include <TBranch.h>
 #include <TBufferJSON.h>
+#include <TMath.h>
 #include "BParticle.h"
 #include "BEvent.h"
 class Hdr{
 };
 void Blab2::h1d(const char *varname, const char *name, int nbins, double min, double max, int id ){
    TString svar(varname);
-   TString axis[]={"mass (GeV/c2)","momentum (GeV/c)","energy (GeV)","charge","identity","momentum (GeV/c)","momentum (GeV/c)","momentum (GeV/c)","momentum (GeV/c)"};
+   TString axis[]={"mass (GeV/c2)",
+                "momentum (GeV/c)",
+                "energy (GeV)","charge",
+                "identity",
+                "momentum (GeV/c)",
+                "momentum (GeV/c)",
+                "momentum (GeV/c)",
+                "momentum (GeV/c)",
+                "angle (deg.)",
+                "cos(theta)"};
    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;
+   if (svar.Contains("GetTheta"))             fHtype[id]=9;
+   if (svar.Contains("GetCosTheta"))          fHtype[id]=10;
    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);
+}
 int Blab2::Fill(int hid, BParticle *p){
-if (hid>=0 && fH[hid]) {
+  if (hid>=0 && fH[hid]) {
-      float val;
+      double val;
       switch (fHtype[hid]){
-      case 0: val = p->GetMass(); break;
+      case 0 : val  = p->GetMass(); break;
-      case 1: val = p->GetMomentum(); break;
+      case 1 : val  = p->GetMomentum(); break;
-      case 2: val = p->e(); break;
+      case 2 : val  = p->e(); break;
-      case 3: val = p->charge(); break;
+      case 3 : val  = p->charge(); break;
-      case 4: val = p->pid(); break;
+      case 4 : val  = p->pid(); break;
-      case 5: val = p->px(); break;
+      case 5 : val  = p->px(); break;
-      case 6: val = p->py(); break;
+      case 6 : val  = p->py(); break;
-      case 7: val = p->pz(); break;
+      case 7 : val  = p->pz(); break;
+      case 8 : val  = p->GetTransverseMomentum(); break;
+      case 9 : val  = (p->GetMomentum()!=0) ? p->pz()/p->GetMomentum() : 0; val = 180*acos(val)/TMath::Pi(); break;
-      case 8: val = sqrt(p->px()*p->px()+p->py()*p->py()); break;
+      case 10: val  = (p->GetMomentum()!=0) ? p->pz()/p->GetMomentum() : 0; break;
-      default : val = 0 ; break;
+      default: val  = 0 ; break;
-      }
+   }
    fH[hid]->Fill(val);
+}
 return 0;
+}
+      }
+   }
  }
  return _p2;
+}
 int Blab2::selector(int pin, int charge, SIMPLEPID type ,  int hid, int pout ){
    hdr.len= msg.Length();
    hdr.progress= progress;
    fwrite(&hdr,3*sizeof(int),1,stdout);
    fwrite(msg.Data(),hdr.len,1,stdout);
    fflush(stdout);
 }
 void Blab2::Process(){
 char sList[0xFFFF];
 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);
 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 i    =TMath::Min(fNfirst, (UInt_t) t-> GetEntries());
 int cNeve=TMath::Min(fNfirst+fNeve, (UInt_t) t-> GetEntries());
+int fPart = fPrint;
+int totaltracks = 0;
 while (i<cNeve){
  t-> GetEntry(i); // Read the content of the event from the file
  fList[0]= mevent->GetParticleList();
  event();
  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++;
    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 (fPart) fHnprt->Fill(nprt);
+ totaltracks += nprt;
  if (fPrint) {
    sprintf(sList,"%s</table>",sList);
    send_message(0, TString(sList),progress);
    nev++;
+ }
  mevent-> Clear();  // Clear the memory.
  for (int k=0;k<100;k++) if (fList[k]!=0) fList[k]->Clear();
  i++;
+}
+double avgtracks=(i)?float(totaltracks)/i:0;
+send_message(0, TString::Format("Number of events processed: %d<br/>\nNumber of particles: %d<br/>\nAverage number of particles per event %f<br/>\n", i, totaltracks, avgtracks ),100);
-send_message(2, TString::Format("Number of events processed %d\n", i ),100);
-send_message(1,TBufferJSON::ConvertToJSON(fHnprt),100 );
+if (fPart) send_message(1,TBufferJSON::ConvertToJSON(fHnprt),100 );
 for (int i=0;i<100;i++) if (fH[i]!=0) send_message(1,TBufferJSON::ConvertToJSON(fH[i]),100 );
 TDatime d;

Subversion Repositories f9daq

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