• Main Page
  • Related Pages
  • Namespaces
  • Data Structures
  • Files
  • File List
  • Globals

beamtime/cern-oct12/go4/AnalysisMacros/spadic_beam_Analysis.cxx (r4864/r4846)

Go to the documentation of this file.
00001 #include <iostream>
00002 #include <string>
00003 #include <typeinfo>
00004 #include <iostream>
00005 #include <vector>
00006 
00007 
00008 #include "Riostream.h"
00009 #include "TRint.h"
00010 #include "TROOT.h"
00011 #include "TStyle.h"
00012 #include "Riostream.h"
00013 #include "TRandom2.h"
00014 
00015 #include "TH2.h"
00016 #include "TH1.h"
00017 #include "TF1.h"
00018 #include "TGaxis.h"
00019 #include "TGraph.h"
00020 #include "TGraphErrors.h"
00021 #include "TMultiGraph.h"
00022 #include "TCanvas.h"
00023 #include "TPad.h"
00024 #include "TProfile.h"
00025 #include "TStopwatch.h"
00026 
00027 #include "TFile.h"
00028 #include "TTree.h"
00029 #include "TBranch.h"
00030 #include "TMath.h"
00031 #include "TLegend.h"
00032 #include "TLine.h"
00033 #include "TColor.h"
00034 #include "TText.h"
00035 #include "TKey.h"
00036 #include "TLatex.h"
00037 #include "TPaveText.h"
00038 #include "TCanvas.h"
00039 
00040 #include "TPrincipal.h"
00041 #include "TMatrix.h"
00042 #include "TFitResultPtr.h"
00043 #include "TFitResult.h"
00044 
00045 #include "TRocEvent.h"
00046 #include "roc/Message.h"
00047 #include "TSpadicEvent.h"
00048 #include "TEpicsEvent.h"
00049 #include "TFiberHodEvent.h"
00050 //#include "TTriglogEvent.h"
00051 #include "TMbsCrateEvent.h"
00052 #include "TBeamMonitorEvent.h"
00053 #include "TCernOct12UnpackEvent.h"
00054 #include "TCernOct12DetectorEvent.h"
00055 #include "spadic_noise_Analysis.cxx"
00056 
00057 void spadic_beam_Analysis(TString filename="data2012/TTrees/Be_run2310004.root", Bool_t debug=false){
00058   TString name, title, outpic;
00059   outpic = filename;
00060   outpic.ReplaceAll("data2012/TTrees/merged/","data2012/pics/EPICS/");
00061   outpic.ReplaceAll(".root","");
00062   
00063   gROOT->Reset(); 
00064   gROOT->SetStyle("Plain");
00065   gStyle->SetPalette(1,0);  
00066   gStyle->SetOptTitle(kFALSE);
00067   gStyle->SetPadTickX(1);                         //<-- tic marks on all axes
00068   gStyle->SetPadTickY(1);  
00069   gStyle->SetFillStyle(4000);
00070   gStyle->SetFrameFillStyle(0);
00071   if (debug){
00072     gStyle->SetOptStat(kFALSE);
00073     gStyle->SetStatStyle(0);
00074   }
00075   gStyle->SetTitleStyle(0);
00076   gStyle->SetCanvasBorderSize(0);
00077   gStyle->SetFrameBorderSize(0);
00078   gStyle->SetLegendBorderSize(0);
00079   gStyle->SetStatBorderSize(0);
00080   gStyle->SetTitleBorderSize(0);
00081   TGraphErrors* gasFlow = new TGraphErrors();
00082   TGraphErrors* temperatur = new TGraphErrors();
00083   TGraphErrors* HVA1_I  = new TGraphErrors();
00084   HVA1_I->SetLineColor(2);
00085   TGraphErrors* HVA1_U = new TGraphErrors();
00086   TGraphErrors* HVD1_I  = new TGraphErrors();
00087   HVD1_I->SetLineColor(2);
00088   TGraphErrors* HVD1_U = new TGraphErrors();
00089   TGraphErrors* HVA2_I = new TGraphErrors();
00090   HVA2_I->SetLineColor(2);
00091   TGraphErrors* HVA2_U = new TGraphErrors();
00092   TGraphErrors* HVD2_I = new TGraphErrors();
00093   HVD2_I->SetLineColor(2);
00094   TGraphErrors* HVD2_U = new TGraphErrors();
00095   TGraphErrors* LVSu1_I = new TGraphErrors();
00096   LVSu1_I->SetLineColor(2);
00097   TGraphErrors* LVSu1_U = new TGraphErrors();
00098   TGraphErrors* LVSp1_I = new TGraphErrors();
00099   LVSp1_I->SetLineColor(2);
00100   TGraphErrors* LVSp1_U = new TGraphErrors();
00101   TGraphErrors* LVSu2_I = new TGraphErrors();
00102   LVSu2_I->SetLineColor(2);
00103   TGraphErrors* LVSu2_U = new TGraphErrors();
00104   TGraphErrors* LVSp2_I = new TGraphErrors();
00105   LVSp2_I->SetLineColor(2);
00106   TGraphErrors* LVSp2_U = new TGraphErrors();
00107  
00108   TFile *inputFile = new TFile(filename,"READ");
00109   if (!inputFile->IsOpen()){
00110     cout << "file not found: " << filename << endl;
00111     return;
00112   } else {
00113     TStopwatch timer;
00114     timer.Start();
00115     TTree* theTree = NULL;
00116     TKey* kee = NULL;
00117     TIter iter(inputFile->GetListOfKeys());
00118     while ( ( kee = dynamic_cast<TKey*>(iter()) ) !=0 ) {
00119       theTree = dynamic_cast<TTree*>(kee->ReadObj());
00120       if (theTree)
00121         break; // we take first Tree in file, disregarding its name...
00122     }
00123     if(theTree == NULL) {
00124       cout <<"Error: no Tree in file "<< filename.Data() << endl;
00125       return;
00126     }
00127     TCernOct12UnpackEvent* evnt = new TCernOct12UnpackEvent;
00128     TGo4EventElement* theBase=evnt;
00129     evnt->synchronizeWithTree(theTree, &theBase);
00130     const Int_t entries = (Int_t)theTree->GetEntries();
00131     printf("%7i Events found in TTree\n",entries);
00132     //TMbsCrateEvent* fCrateInputEvent = NULL;
00133     //TSpadicEvent* SpadicInputEvent = NULL;
00134     //TSpadicData* theSpadic1 = NULL;
00135     //TSpadicData* theSpadic2 = NULL;
00136     TEpicsEvent* EpicsInputEvent = NULL;
00137     Int_t j = 0;
00138     for(Int_t i = 0; i < entries; ++i) {
00139       Statusbar(i,entries);
00140       theTree->GetEntry(i);    
00141       //fCrateInputEvent = dynamic_cast<TMbsCrateEvent*>(evnt->GetSubEvent("MBSCRATE"));    
00142       //SpadicInputEvent = dynamic_cast<TSpadicEvent*>(evnt->GetSubEvent("SPADIC"));
00143       //theSpadic1 = dynamic_cast<TSpadicData*>(SpadicInputEvent->getEventElement(usedSuId[0]));//0
00144       //theSpadic2 = dynamic_cast<TSpadicData*>(SpadicInputEvent->getEventElement(usedSuId[1]));//1  
00145       EpicsInputEvent=dynamic_cast<TEpicsEvent*>(evnt->GetSubEvent("EPICS"));
00146       //Ch1 = fCrateInputEvent->fData1182[1][0];
00147       //Ch2 = fCrateInputEvent->fData1182[0][0];
00148       //Pb  = fCrateInputEvent->fData1182[1][1];
00149       //gasFlow.push_back(EpicsInputEvent->GetDouble("CBM:TRD:GFM0:Flow"));
00150      
00151       if (EpicsInputEvent->GetDouble("CBM:TRD:GFM0:Flow") > 0){
00152         j = gasFlow->GetN();
00153         gasFlow->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:GFM0:Flow"));
00154         //j = LVSu1_U->GetN();
00155         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U001:vmon") > 0)
00156         LVSu1_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U001:vmon"));
00157         //j = LVSu1_I->GetN();
00158         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U001:imon") > 0)
00159         LVSu1_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U001:imon"));
00160         //j = LVSp1_U->GetN();
00161         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U000:vmon") > 0)
00162         LVSp1_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U000:vmon"));
00163         //j = LVSp1_I->GetN();
00164         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U000:imon") > 0)
00165         LVSp1_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U000:imon"));
00166         //j = HVA1_U->GetN();
00167         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U300:vmon") > 0)
00168         HVA1_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U300:vmon"));
00169         //j = HVA1_I->GetN();
00170         // if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U300:imon") > 0)
00171         HVA1_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U300:imon") * 1e9);
00172         //j = HVD1_U->GetN();
00173         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U200:vmon") > 0)
00174         HVD1_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U200:vmon"));
00175         //j = HVD1_I->GetN();
00176         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U200:imon") > 0)
00177         HVD1_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U200:imon") * 1e9);
00178 
00179         //j = LVSu2_U->GetN();
00180         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U003:vmon") > 0)
00181         LVSu2_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U003:vmon"));
00182         //j = LVSu2_I->GetN();
00183         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U003:imon") > 0)
00184         LVSu2_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U003:imon"));
00185         //j = LVSp2_U->GetN();
00186         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U002:vmon") > 0)
00187         LVSp2_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U002:vmon"));
00188         //j = LVSp2_I->GetN();
00189         //if (EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U002:imon") > 0)
00190         LVSp2_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:WIENLV_U002:imon"));
00191         //j = HVA2_U->GetN();
00192         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U301:vmon") > 0)
00193         HVA2_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U301:vmon"));
00194         //j = HVA2_I->GetN();
00195         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U301:imon") > 0)
00196         HVA2_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U301:imon") * 1e9);
00197         //j = HVD2_U->GetN();
00198         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U201:vmon") > 0)
00199         HVD2_U->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U201:vmon"));
00200         //j = HVD2_I->GetN();
00201         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U201:imon") > 0)
00202         HVD2_I->SetPoint(j,Double_t(i)/1000,EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U201:imon") * 1e9);
00203         //j = temperatur->GetN();
00204         //if (EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U900:temp") > 0)
00205         temperatur->SetPoint(j,Double_t(i)/1000,(
00206                                                  //EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U000:temp")+
00207                                                  //EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U100:temp")+
00208                                                  //EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U200:temp")+
00209                                                  //EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U300:temp")+
00210                                                  //EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U400:temp")+
00211                                                  EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U500:temp")+
00212                                                  EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U600:temp")+
00213                                                  EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U700:temp")+
00214                                                  EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U800:temp")+
00215                                                  EpicsInputEvent->GetDouble("CBM:TRD:ISEGHV_U900:temp")) / 5.
00216                              );
00217       }
00218     }
00219     TMultiGraph *HVA1 = new TMultiGraph();
00220     HVA1->Add(HVA1_U);
00221     HVA1->Add(HVA1_I);
00222     TMultiGraph *HVD1 = new TMultiGraph();
00223     HVD1->Add(HVD1_U);
00224     HVD1->Add(HVD1_I);
00225     TMultiGraph *HVA2 = new TMultiGraph();
00226     HVA2->Add(HVA2_U);
00227     HVA2->Add(HVA2_I);
00228     TMultiGraph *HVD2 = new TMultiGraph();
00229     HVD2->Add(HVD2_U);
00230     HVD2->Add(HVD2_I);
00231     TMultiGraph *LVSu1 = new TMultiGraph();
00232     TMultiGraph *LVSu2 = new TMultiGraph();
00233     TMultiGraph *LVSp1 = new TMultiGraph();
00234     TMultiGraph *LVSp2 = new TMultiGraph();
00235     TCanvas *c = new TCanvas("c","c",2*800,600);
00236     TCanvas *c1 = new TCanvas("c1","c1",2*800,2*600);
00237     c1->Divide(2,2);
00238     c->Divide(2,1);
00239     c->cd(1);
00240     gasFlow->Draw("AL");
00241     gasFlow->SetMinimum(5);
00242     gasFlow->SetMaximum(13);
00243     gasFlow->GetXaxis()->SetTitle("event ID/1000");
00244     gasFlow->GetYaxis()->SetTitle("gas flow [l/h]");
00245     
00246     c->cd(2);
00247     temperatur->Draw("AL");
00248     temperatur->SetMinimum(15);
00249     temperatur->SetMaximum(30);
00250     temperatur->GetXaxis()->SetTitle("event ID/1000");
00251     temperatur->GetYaxis()->SetTitle("temperature [{}^{o}C]");
00252     c1->cd(1);
00253     //HVA1_U->Draw("AL");
00254     HVA1->Draw("AL");
00255     HVA1->SetMinimum(0);
00256     HVA1->SetMaximum(2000);
00257     HVA1->GetXaxis()->SetTitle("event ID/1000");
00258     HVA1->GetYaxis()->SetTitle("anode voltage [V]");
00259     HVA1->GetYaxis()->SetTitleOffset(1.3);
00260     TGaxis *axisA = new TGaxis(c1->cd(1)->GetUxmax(),c1->cd(1)->GetUymin(),
00261                                c1->cd(1)->GetUxmax(),c1->cd(1)->GetUymax(),
00262                                HVA1->GetYaxis()->GetXmin(),HVA1->GetYaxis()->GetXmax(),
00263                                510,"+L");
00264     axisA->SetLineColor(kRed);
00265     axisA->SetTextColor(kRed);
00266     axisA->SetLabelColor(kRed);
00267     axisA->SetTitleOffset(1.3);
00268     axisA->SetTitle("anode current [nA]");
00269     axisA->Draw("");
00270     c1->cd(2);
00271     //HVA1_I->Draw("AL");
00272     HVD1->Draw("AL");
00273     HVD1->SetMinimum(0);
00274     HVD1->SetMaximum(750);
00275     HVD1->GetXaxis()->SetTitle("event ID/1000");
00276     HVD1->GetYaxis()->SetTitle("drift voltage [V]");
00277     HVD1->GetYaxis()->SetTitleOffset(1.3);
00278     TGaxis *axisD = new TGaxis(c1->cd(2)->GetUxmax(),c1->cd(2)->GetUymin(),
00279                                c1->cd(2)->GetUxmax(),c1->cd(2)->GetUymax(),
00280                                HVD1->GetYaxis()->GetXmin(),HVD1->GetYaxis()->GetXmax(),
00281                                510,"+L");
00282     axisD->SetLineColor(kRed);
00283     axisD->SetTextColor(kRed);
00284     axisD->SetLabelColor(kRed);
00285     axisD->SetTitleOffset(1.3);
00286     axisD->SetTitle("drift current [nA]");
00287     axisD->Draw("");
00288     c1->cd(3);
00289     //HVD1_U->Draw("AL");
00290     HVA2->Draw("AL");
00291     HVA2->SetMinimum(0);
00292     HVA2->SetMaximum(2000);
00293     HVA2->GetXaxis()->SetTitle("event ID/1000");
00294     HVA2->GetYaxis()->SetTitle("voltage [V]");
00295     HVA2->GetYaxis()->SetTitleOffset(1.3);
00296     axisA->Draw("");
00297     c1->cd(4);
00298     //HVD1_I->Draw("AL");
00299     HVD2->Draw("AL");
00300     HVD2->SetMinimum(0);
00301     HVD2->SetMaximum(750);
00302     HVD2->GetXaxis()->SetTitle("event ID/1000");
00303     HVD2->GetYaxis()->SetTitle("voltage [V]");
00304     HVD2->GetYaxis()->SetTitleOffset(1.3);
00305     axisD->Draw("");
00306     c->SaveAs(outpic+".pdf");
00307     c1->SaveAs(outpic+"_HV.pdf");
00308     c->SaveAs(outpic+".png");
00309     c1->SaveAs(outpic+"_HV.png");
00310     /*
00311       c1->cd(5);
00312       LVSu1_U->Draw("AL");
00313       LVSu1_U->SetMinimum(0);
00314       LVSu1_U->SetMaximum(10);
00315       c1->cd(6);
00316       LVSu1_I->Draw("AL");
00317       LVSu1_I->SetMinimum(0);
00318       LVSu1_I->SetMaximum(10);
00319       c1->cd(7);
00320       LVSp1_U->Draw("AL");
00321       LVSp1_U->SetMinimum(0);
00322       LVSp1_U->SetMaximum(10);
00323       c1->cd(8);
00324       LVSp1_I->Draw("AL");
00325       LVSp1_I->SetMinimum(0);
00326       LVSp1_I->SetMaximum(10);
00327       c1->cd(9);
00328       //HVA2_U->Draw("AL");
00329       c1->cd(10);
00330       //HVA2_I->Draw("AL");
00331       c1->cd(11);
00332       //HVD2_U->Draw("AL");
00333       c1->cd(12);
00334       //HVD2_I->Draw("AL");
00335       c1->cd(13);
00336       LVSu2_U->Draw("AL");
00337       c1->cd(14);
00338       LVSu2_I->Draw("AL");
00339       c1->cd(15);
00340       LVSp2_U->Draw("AL");
00341       c1->cd(16);
00342       LVSp2_I->Draw("AL");
00343     */
00344     timer.Stop();
00345     Double_t rtime = timer.RealTime();
00346     Double_t ctime = timer.CpuTime();
00347     printf("\n\n*********************************************************\n");
00348     //printf("   RealTime=%f seconds,  CpuTime=%f seconds\n",rtime,ctime);
00349     printf("    RealTime = %i min. %.1f s,  CpuTime = %i min. %.1f s\n",Int_t(rtime/60),(rtime/60 - Int_t(rtime/60))*60.,Int_t(ctime/60),(ctime/60 - Int_t(ctime/60))*60.);
00350     printf("             %.1f event/s,             %.1f event/s\n",Float_t(entries)/rtime,Float_t(entries)/ctime);
00351     printf("*********************************************************\n\n");
00352     inputFile->Close();
00353   }
00354 }

Generated on Tue Dec 10 2013 04:52:18 for ROCsoft by  doxygen 1.7.1