/[PAMELA software]/calo/flight/CaloBragg/doc/example.cpp
ViewVC logotype

Diff of /calo/flight/CaloBragg/doc/example.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1 by pamelats, Fri Jun 13 08:11:03 2008 UTC revision 1.2 by mocchiut, Thu Aug 19 11:55:32 2010 UTC
# Line 1  Line 1 
 //Nella release 4 ci sono i TGraph e i fit sui TGraph. Ho perso così l'info in Z che avevo.  
   
   
1  #include <fstream>  #include <fstream>
2  #include <sstream>  #include <sstream>
3  #include <string>  #include <string>
# Line 33  Line 30 
30  #include <TCanvas.h>  #include <TCanvas.h>
31  #include <PamLevel2.h>  #include <PamLevel2.h>
32    
33  #include <CaloNuclei.h> // ADD HEADER FILE HERE  #include <CaloBragg.h> // ADD HEADER FILE HERE
34    
35  using namespace std;  using namespace std;
36    
37    
38    //===================================
39  void example(TString dir, TString flist, TString treelist="+ALL"){  // global variables
40    //===================================
41    Bool_t    DEBUG;
42    Bool_t    SIMU;
43    TString   DIR;
44    TString   LIST;
45    TString   OPTIONS;
46    ULong64_t MAXEV;
47    TString   OUTFILE;
48    
49    PamLevel2 *pam_event = NULL;
50    TFile     *outfv = NULL;
51    
52    bool fillTree   = false;
53    bool fillHistos = false;
54    
55    
56    //==========================================
57    //000000000000000000000000000000000000000000
58    //==========================================
59    Int_t Loop(TString ddir,TString list, ULong64_t nmax, TString options, TString outfile){
60    
61        
62    PamLevel2 *pam_event = new PamLevel2(dir,flist);       // << create pamela event    TH1F *hzeta = new TH1F("hzeta","Zeta Distribution",100,0.,50.);
63    TChain *T = pam_event->GetPamTree();    
64      if(options.Contains("fillTree"))fillTree=true;
65    CaloNuclei *cn = new CaloNuclei(pam_event);  // CREATE CaloNuclei OBJECT GIVING PamLevel2 POINTER    if(options.Contains("fillHisto"))fillHistos=true;
66      
67      //  -------------------
68    Int_t nevent = pam_event->GetEntries();    //  create output files
69    cout<< " Events= "<<nevent<<"\n";    //  -------------------
70      TString outfile_t ="";
71      //  TString outfile_h ="";
72      //  TString outfile_v ="";
73      if( outfile.IsNull() ){
74        if(!list.IsNull())outfile = list(0,list.Last('.'));    
75        else outfile = "output";        
76      }else{        
77        if(outfile.Contains(".root"))outfile = outfile(0,outfile.Last('.'));
78      }
79      //    if(fillTree){
80      outfile_t = outfile;
81      outfile_t.Append("-hzeta.root");
82      outfv = (TFile*)gROOT->FindObject(outfile_t);
83      if (outfv) outfv->Close();
84      outfv = new TFile(outfile_t,"RECREATE");
85      if(outfv->IsZombie()){
86        cout << "Output file could not be created\n";
87        return 1;
88      };
89      cout << "Created output file: "<<outfv->GetName()<<endl;
90      //  --------------------
91      //  read input file/list
92      //  --------------------
93      pam_event = new PamLevel2(ddir,list,options);
94      //
95      CaloBragg *cn = new CaloBragg(pam_event);  // CREATE CaloBragg OBJECT GIVING PamLevel2 POINTER
96    
   const Int_t size=nevent;  
   
97    //    //
98    // loop over events    // loop over events
99    //    //
100    //  for (Int_t i=0; i<nevent;i++){    ULong64_t nevents = pam_event->GetEntries();  
101    for (Int_t i=0; i<50;i++){    if(!nmax)nmax = numeric_limits<ULong64_t>::max();    
102          if(nevents < nmax)nmax=nevents;  
103      cout << endl<<" Start loop over events:  "<< nmax<<endl;
104      //    
105      for(ULong64_t i=0; i<nmax; i++){
106        //
107      if ( i%10000 == 0) cout<<i<<endl;      if ( i%10000 == 0) cout<<i<<endl;
108      //      //
109      T->GetEntry(i);  // get the entry, also pam_event->GetEntry is fine      pam_event->Clear();
110            if( pam_event->GetEntry(i) ){ // get the entry
     cn->Print(); // this will simply print on the terminal the determined CaloNuclei values  
       
     Float_t mip1 = cn->Get_dEdx1();    // here how you can use the "old" dE/dx on the first plane  
       
     if ( cn->Get_qpremeanN() > 0. && cn->Get_interplane()>9 ) { // other new variables  
111        //        //
112        // do something        cn->Print(); // this will simply print on the terminal the determined CaloNuclei values
113        //        //    
114      };        Float_t zeta = cn->Get_lpz();    // here how you can use the "old" dE/dx on the first plane
115          printf(" Event number %llu zeta %f \n",i,zeta);
116          hzeta->Fill(zeta);
117    
118        }  
119      };
120      pam_event->Clear();
121      //
122      // save histo to file
123      //
124      gROOT->cd();
125      outfv->cd();
126      //
127      hzeta->Write();
128      outfv->Close();
129      //
130      printf(" end!\n");
131      return 0;
132        
133    }
134    
135      // if ( nrig > 5. ) caloz->Fill(mip1);  /////////////////////////////////////////////////////////////////
136    
137      //  void usage(){
     // if you want to do it for other tracks and not only for the track number 0 (default) you can do this way:  
     //  
     Float_t mipt[100];  
     memset(mipt,0.,100*sizeof(Float_t));  
     for (Int_t nt = 0; nt < pam_event->GetTrkLevel2()->GetNTracks(); nt++){  
       PamTrack *track = pam_event->GetTrack(nt);  
       cn->Process(nt);  
       if ( nt < 100 ) mipt[nt] = cn->Get_dEdx1();    
       if ( nt > 0 ){  
         cn->Print();  
       };  
     };  
     //  
138    
139    }      cout << "------------------------------------------------------------"<<endl;
140    printf(" end!\n");    cout << "Loop over events (on one or more Level2-files), applying some selection cuts (defined in My-Selection.cpp), \n";
141      cout << "creating output histograms (defined in My-Histos.cpp) and/or trees with selected events. \n \n ";
142  };    cout << "USAGE:"<<endl;
143      cout << "-processDir  DIR     -  Level2 data directory \n";
144      cout << "-processList LIST    -  list  of files (.txt) or single file (.root) to be analysed \n";
145      cout << "-outputFile  PATH    -  name of the output file \n";
146      cout << "-NumEvents   XXX     -  number of events to be analysed \n";
147      cout << "--debug, -g          -  debug mode \n";
148      cout << "--help, -h           -  print this help \n";
149      cout << "-options [ options ] -  options: \n";
150      cout << "                        fillHistos --> create an output file with histograms  \n";
151      cout << "                        fillTree   --> create an output file with trees storing the selected events \n ";
152      cout << "                        +(-)ALL    --> inlcude(exclude) all trees and branches \n "   ;
153      cout << "                        +(-)TRK1 +(-)TRK2 +(-)CAL1 +(-)CAL2 +(-)TOF +(-)TRG +(-)ND +(-)S4 +(-)ORB --> inlcude(exclude) trees and branches  \n"  ;
154      cout << "------------------------------------------------------------"<<endl;
155    }
156  //  //
157    int HandleInputPar(int argc, char **argv){
158    
159      if(argc>1){
160                    
161        if(!strcmp(argv[1], "-h") || !strcmp(argv[1], "--help") ){
162          usage();
163          return(1);
164        };
165        // -----------------------
166        // Read input parameters
167        // -----------------------  
168        DEBUG   = false;
169        SIMU   = false;
170        DIR     = gSystem->WorkingDirectory();
171        LIST    = "";
172        OUTFILE = "";
173        OPTIONS = "+AUTO fillTree";
174        MAXEV   = 0;        
175                    
176    
177        for (int i = 1; i < argc; i++){            
178          // -----------------------------------------------------//
179          if (!strcmp(argv[i], "-processDir")){
180            if (++i >= argc) throw -1;
181            DIR = argv[i];
182            cout << "processDir "<<DIR<<endl;
183            continue;
184          }
185          // -----------------------------------------------------//
186          else if (!strcmp(argv[i], "-processList")){
187            if (++i >= argc) throw -1;
188            LIST = argv[i];
189            cout << "processList "<<LIST<<endl;
190            continue;
191          }  
192          // -----------------------------------------------------//
193          else if (!strcmp(argv[i], "-outputFile")){
194            if (++i >= argc) throw -1;
195            OUTFILE = argv[i];
196            cout << "outputFile "<<OUTFILE<<endl;
197            continue;
198          }  
199          // -----------------------------------------------------//
200          else if (!strcmp(argv[i], "-NumEvents")){
201            if (++i >= argc) throw -1;
202            MAXEV = atoi(argv[i]);
203            cout << "NumEvents "<<MAXEV<<endl;
204            continue;      
205          }
206          // -----------------------------------------------------//
207          else if (!strcmp(argv[i], "-options")){
208            if (++i >= argc) throw -1;
209            OPTIONS = argv[i];
210            if( OPTIONS.Contains("[") ){
211              do{
212                if (++i >= argc) throw -1;
213                OPTIONS.Append(argv[i]);
214              }while(!OPTIONS.Contains("]"));
215            }else cout << "wrong option format --> ignoring " << endl;
216          }
217          // -----------------------------------------------------//
218          else if (!strcmp(argv[i], "--debug") || !strcmp(argv[i], "-g")){
219            DEBUG = true;
220            continue;      
221          }
222          else if (!strcmp(argv[i], "--simu") || !strcmp(argv[i], "-s")){
223            SIMU = true;
224            continue;      
225          }
226          // -----------------------------------------------------//
227          else{
228            cout << "Unidentified input parameter. Ignored."<< endl;
229          };
230        };
231      }else{
232        usage();
233        return(1);
234      };
235      // -----------------------
236      // Check input parameters
237      // -----------------------  
238            
239                    
240      return(0);
241            
242    };
243    //
244    
245    int main(int argc, char **argv)
246    {
247            
248      if( HandleInputPar(argc,argv) )return(1);
249            
250      //    Loop(DIR,LIST,MAXEV,"-ALL+TRK1+TRK2+CAL1+CAL2+TOF+AC",OUTFILE);
251      cout << "OPTIONS "<<OPTIONS<<endl;
252      Int_t tt = Loop(DIR,LIST,MAXEV,OPTIONS,OUTFILE);
253        
254      cout << "Back to main - end"<<endl;
255        
256      return tt;
257            
258    };
259    
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.2
  ViewVC Help
Powered by ViewVC 1.1.23