/[PAMELA software]/quicklook/QLflightTmtc_Header/HeaderScan.cpp
ViewVC logotype

Annotation of /quicklook/QLflightTmtc_Header/HeaderScan.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.1.1.2 - (hide annotations) (download) (vendor branch)
Sat Jun 17 11:02:09 2006 UTC (18 years, 5 months ago) by pam-rm2
Branch: QLflightTmtc_Header
CVS Tags: second, v1r04, v1r03, v1r1
Changes since 1.1.1.1: +217 -206 lines
aggiornati anche il PacketScan e l'HeaderScan

1 pam-rm2 1.1 /**
2     * Header Scan
3     * Author Nagni
4     * version 1.0
5     *
6     * Version 1.1 - 28 December 2004
7     * If outList does not exist the function exit to prompt
8     * If nevents = 0 the function exit to promt
9     *
10     * Version 1.2 - 3 January 2005
11     * Two canvases are created to see the graphs better
12     *
13     * Version 1.25 - 13 January 2005
14     * Changed Int_t to Float because of variable range size
15     * (UInt_t has been excluded beacuse of uncompatibility with TGraph)
16     *
17     * version 1.3 - 22 February 2005 - Nagni
18     * For compatibility with batch mode excution:
19     * 1) Added "include <iostream>" and "using namespace std"
20     * 2) Removed gROOT->Reset()
21     *
22     * Version 1.4
23     * Date 08 March 2005
24     * Added "format" parameter for saving the produced image in various formats
25     * (for a complete list of types refer to TPad::SaveAs method)
26     *
27     * Version 1.5
28     * Date 09 February 2006 - Marcelli
29     * Update to work with new Yoda output
30     *
31     * Version 1.6
32     * Date 27 February 2006 - Marcelli
33     * For compilation:
34     * Added function "int main(int argc, char* argv[])"
35     *
36     *
37     * Description: This script creates two canvases with five pads. The first pad shows packetID variable (for all packets) vs. OBT.
38     * The second pad shows the number of physic packets vs. OBT. The third pad shows the lenght of Physic packets (byte) vs. OBT.
39     * The fourth pad shows the packetcounter of physic packets vs. OBT. The fifth pad shows PacketCounter vs. File Offset.
40     *
41     * Parameters:
42     * TSTring base - the path to the root directory for the specific Pamela unpack session
43     * There is no default value, without this input the program will not run
44     * TString outDir - the path where to save the output image (Default = ./)
45     * TString format - the format which will be used fo rsave the produced images (Default = "jpg")
46 pam-rm2 1.1.1.2 *
47     *
48     * Version 1.7
49     * Date 16 June 2006 - Malvezzi
50     *
51     * Description of changes:
52     * Implementation of case: numebr of events are <= 0.
53     * Remove of the graph "grPcktId1"; see PacketScan for the same information.
54     * Fixed bugs: for a large namber of events is not possible to have vectors, so all graphs have been converted in histograms
55     *
56 pam-rm2 1.1 */
57    
58    
59 pam-rm2 1.1.1.2 #include <fstream>
60     #include <math.h>
61     #include "TLatex.h"
62     #include "TF1.h"
63     #include "TPaveText.h"
64     #include "TMultiGraph.h"
65 pam-rm2 1.1 #include <sstream>
66     #include <iostream>
67     #include "TString.h"
68     #include "TStyle.h"
69     #include "TFile.h"
70     #include "TList.h"
71     #include "TTree.h"
72     #include "TObjString.h"
73     #include "TCanvas.h"
74     #include "TGraph.h"
75     #include "TH1F.h"
76     #include "TGaxis.h"
77     #include "EventHeader.h"
78     #include "PscuHeader.h"
79    
80     using namespace std;
81    
82     void HeaderScan(TString base, TString outDir, TString format)
83     {
84    
85 pam-rm2 1.1.1.2 //------------------- Variables initilization -------------------------//
86    
87 pam-rm2 1.1 TList *list = new TList;
88     Int_t numkey;
89     TObject *key = new TObject;
90     const char *name;
91 pam-rm2 1.1.1.2 //TTree* tr = new TTree;
92     //TTree* tr1 = new TTree;
93     Long64_t nevents=0; // ev=0, events=0;
94     ULong_t lastime, firstime;
95     double obmin=0.;
96     double obmax=0.;
97     double obt;
98     ULong_t maxoffset, minoffset;
99     Float_t pcktLenght =0., pcktCounter=0., offset=0., allCounter=0.;
100     stringstream oss, oss1, oss2, oss3, noentries;
101 pam-rm2 1.1 //------- load root file --------------
102     TFile *file = new TFile(base.Data());
103     if (!file){
104     printf("No such file in the directory has been found");
105     return;
106     }
107     if (outDir == "" ) outDir = ".";
108 pam-rm2 1.1.1.2 list = file->GetListOfKeys(); //get list of trees in the file
109     numkey = file->GetNkeys(); //get number of trees in the file
110 pam-rm2 1.1
111 pam-rm2 1.1.1.2 TTree *PhysicsTr = (TTree*)file->Get("Physics");
112     TBranch *headBr = PhysicsTr->GetBranch("Header");
113    
114     pamela::EventHeader *eh = 0;
115     pamela::PscuHeader *ph = 0;
116    
117     PhysicsTr->SetBranchAddress("Header", &eh);
118    
119     nevents = headBr->GetEntries();
120     const Int_t sizetot = nevents;
121    
122     TString filename = ((TObjString*)base.Tokenize('/')->Last())->GetString();
123     filename = ((TObjString*)filename.Tokenize('.')->First())->GetString();
124    
125     //----------- If nevents < = 0 ---------------------------------/
126     if (nevents<=0) {
127     printf("nevents = %i \n", nevents);
128     printf(" \n");
129    
130     TCanvas *canv = new TCanvas("No entries", "No entries ", 400, 200);
131     canv->SetFillColor(10);
132     canv->cd();
133    
134     TLatex *l = new TLatex();
135     l->SetTextAlign(12);
136     l->SetTextSize(0.15);
137     l->SetTextColor(2);
138     noentries.str("");
139     noentries<< "HeaderScan_QL:";
140     l->DrawLatex(0.05, 0.7, noentries.str().c_str());
141     noentries.str("");
142     noentries<< "No Physics entries for this files";
143     l->DrawLatex(0.05, 0.5, noentries.str().c_str());
144    
145     if (outDir == "./") {
146     oss.str("");
147     oss << filename.Data() << "_HeaderScan_QL." << format.Data();
148     } else {
149     oss.str("");
150     oss << outDir.Data() << filename.Data() << "_HeaderScan_QL." << format.Data();
151     }
152 pam-rm2 1.1
153 pam-rm2 1.1.1.2 canv->Update();
154     canv->SaveAs(oss.str().c_str());
155    
156     return;
157     }
158    
159    
160     ///-------------- to know the max and min File Offset ----------------------------//
161     /*for (Int_t i=0; i<numkey; i++){
162 pam-rm2 1.1 key = list->At(i);
163     name=(key->GetName());
164     tr = (TTree*)file->Get(name);
165     if (tr->IsZombie()) continue;
166 pam-rm2 1.1.1.2 events = tr->GetEntries();
167 pam-rm2 1.1 tr->SetBranchAddress("Header", &eh);
168 pam-rm2 1.1.1.2 for (Int_t j = 0; j < events; j++){
169 pam-rm2 1.1 tr->GetEntry(j);
170     ph = eh->GetPscuHeader();
171 pam-rm2 1.1.1.2 ph->GetFileOffset();
172     if(ph->GetFileOffset() <= minoffset) minoffset=ph->GetFileOffset();
173     if(ph->GetFileOffset() >= maxoffset) maxoffset=ph->GetFileOffset();
174     }
175     }*/
176     //*************************** Histograms ************************************************************//
177     //------------------------ First histogram -----------------------------------//
178     headBr->GetEntry(0);
179     ph = eh->GetPscuHeader();
180     firstime = ph->GetOrbitalTime();
181     headBr->GetEntry(nevents-1);
182     ph = eh->GetPscuHeader();
183     lastime = ph->GetOrbitalTime();
184     obmin=firstime;
185     obmax=lastime;
186     oss1.str("");
187     oss1 << filename.Data() << ": Physics Packet per minute;" <<" start @ " << firstime << ", end @ "<< lastime <<"ms";
188     Int_t nbin = (lastime-firstime)/60000;
189     TH1F *h1 = new TH1F ("histo1", oss1.str().c_str(), nbin, obmin, obmax);
190    
191     //------------------------ Second histogram -----------------------------------//
192     oss2.str("");
193     oss2 << filename.Data() << ": Lenght of Physic packets;";
194     Int_t nint = (lastime-firstime);
195     TH1F *PcktLenght = new TH1F ("histo2", oss2.str().c_str(), nint, obmin, obmax);
196    
197     //------------------------ Third histogram -----------------------------------//
198     oss3.str("");
199     oss3 << filename.Data() << ": Physics Counter";
200     TH1F *PcktCounter = new TH1F ("histo3", oss3.str().c_str(), nint, obmin, obmax);
201    
202     //------------------------ Fourth histogram -----------------------------------//
203     /*oss3.str("");
204     oss3 << filename.Data() << ": Packet Counter";
205     Int_t nintoffset = (maxoffset-minoffset)/10;
206     TH1F *FileOffset = new TH1F ("histo4", oss3.str().c_str(), nintoffset, minoffset, maxoffset);
207     */
208     //**************************************************************************************************//
209     //------- fill histograms ---------//
210    
211     for (Int_t i = 0; i < nevents; i++){
212     headBr->GetEntry(i);
213     ph = eh->GetPscuHeader();
214     obt = ph->GetOrbitalTime();
215     pcktLenght=ph->GetPacketLenght();
216     pcktCounter=ph->GetCounter();
217     h1->Fill(obt);
218     PcktLenght->Fill(obt,pcktLenght);
219     PcktCounter->Fill(obt,pcktCounter);
220 pam-rm2 1.1 }
221    
222 pam-rm2 1.1.1.2 /*for (Int_t i=0; i<numkey; i++){
223 pam-rm2 1.1 key = list->At(i);
224     name=(char *)(key->GetName());
225 pam-rm2 1.1.1.2 tr1 = (TTree*)file->Get(name);
226     if (tr1->IsZombie()) continue;
227     tr1->SetBranchAddress("Header", &eh);
228     ev = tr1->GetEntries();
229     for (Int_t j = 0; j < ev; j++){
230     tr1->GetEntry(j);
231 pam-rm2 1.1 ph = eh->GetPscuHeader();
232 pam-rm2 1.1.1.2 offset=ph->GetFileOffset();
233     allCounter= ph->GetCounter();
234     FileOffset->Fill(offset,allCounter);
235 pam-rm2 1.1 }
236 pam-rm2 1.1.1.2 }*/
237    
238     //****************************** Canvases *******************************//
239     //TGaxis::SetMaxDigits(4);
240     //------------------- First Canvas --------------------------------//
241 pam-rm2 1.1 TCanvas *finalCanv1 = new TCanvas("Header_1", base, 1280, 1024);
242 pam-rm2 1.1.1.2 finalCanv1->Divide(1,3);
243 pam-rm2 1.1 finalCanv1->SetFillColor(10);
244    
245     finalCanv1->cd(1);
246 pam-rm2 1.1.1.2 h1->SetStats(kFALSE);
247     h1->GetXaxis()->SetTitle("OBT (ms)");
248 pam-rm2 1.1 h1->GetXaxis()->CenterTitle();
249 pam-rm2 1.1.1.2 h1->GetXaxis()->SetLabelSize(0.04);
250     h1->GetYaxis()->SetTitle("Number of events ");
251 pam-rm2 1.1 h1->GetYaxis()->CenterTitle();
252 pam-rm2 1.1.1.2 h1->GetYaxis()->SetTitleSize(0.06);
253     h1->GetYaxis()->SetTitleOffset(0.8);
254 pam-rm2 1.1 h1->Draw();
255    
256     finalCanv1->cd(2);
257 pam-rm2 1.1.1.2 PcktLenght->SetStats(kFALSE);
258     PcktLenght->GetXaxis()->SetTitle("OBT (ms)");
259     PcktLenght->GetXaxis()->CenterTitle();
260     PcktLenght->GetXaxis()->SetLabelSize(0.04);
261     PcktLenght->GetYaxis()->SetTitle("Lenght (byte)");
262     PcktLenght->GetYaxis()->CenterTitle();
263     PcktLenght->GetYaxis()->SetLabelSize(0.04);
264     PcktLenght->GetYaxis()->SetTitleSize(0.06);
265     PcktLenght->GetYaxis()->SetTitleOffset(0.8);
266     PcktLenght->SetMarkerColor(2);
267     PcktLenght->SetMarkerSize(.5);
268     PcktLenght->SetMarkerStyle(21);
269     PcktLenght->Draw("9p");
270    
271     finalCanv1->cd(3);
272     PcktCounter->SetStats(kFALSE);
273     PcktCounter->GetXaxis()->SetTitle("OBT (ms)");
274     PcktCounter->GetXaxis()->SetTitleSize(0.05);
275     PcktCounter->GetXaxis()->CenterTitle();
276     PcktCounter->GetXaxis()->SetLabelSize(0.04);
277     PcktCounter->GetYaxis()->SetTitle("Counter");
278     PcktCounter->GetYaxis()->SetTitleSize(0.05);
279     PcktCounter->GetYaxis()->CenterTitle();
280     PcktCounter->GetYaxis()->SetLabelSize(0.04);
281     PcktCounter->GetYaxis()->SetTitleSize(0.06);
282     PcktCounter->GetYaxis()->SetTitleOffset(0.8);
283     PcktCounter->SetMarkerColor(4);
284     PcktCounter->SetMarkerSize(.5);
285     PcktCounter->SetMarkerStyle(21);
286     PcktCounter->Draw("9p");
287     //--------- Second Canvas -----------------------------//
288     /*TCanvas *finalCanv2 = new TCanvas("Header_2", base, 1280, 1024);
289     finalCanv2->Divide(1,2);
290     finalCanv2->SetFillColor(10);
291    
292 pam-rm2 1.1 finalCanv2->cd(1);
293 pam-rm2 1.1.1.2 PcktCounter->SetStats(kFALSE);
294     PcktCounter->GetXaxis()->SetTitle("OBT (ms)");
295     PcktCounter->GetXaxis()->SetTitleSize(0.05);
296     PcktCounter->GetXaxis()->CenterTitle();
297     PcktCounter->GetXaxis()->SetLabelSize(0.04);
298     PcktCounter->GetYaxis()->SetTitle("Counter");
299     PcktCounter->GetYaxis()->SetTitleSize(0.05);
300     PcktCounter->GetYaxis()->CenterTitle();
301     PcktCounter->GetYaxis()->SetLabelSize(0.04);
302     PcktCounter->SetMarkerColor(4);
303     PcktCounter->SetMarkerSize(.5);
304     PcktCounter->SetMarkerStyle(21);
305     PcktCounter->Draw("9p");
306 pam-rm2 1.1
307     finalCanv2->cd(2);
308 pam-rm2 1.1.1.2 FileOffset->SetStats(kFALSE);
309     FileOffset->GetXaxis()->SetTitle("File Offset");
310     FileOffset->GetXaxis()->CenterTitle();
311     FileOffset->GetXaxis()->SetTitleSize(0.05);
312     FileOffset->GetXaxis()->SetLabelSize(0.04);
313     FileOffset->GetYaxis()->SetTitle("Counter");
314     FileOffset->GetYaxis()->CenterTitle();
315     FileOffset->GetYaxis()->SetTitleSize(0.05);
316     FileOffset->GetYaxis()->SetLabelSize(0.04);
317     FileOffset->SetMarkerColor(3);
318     FileOffset->SetMarkerSize(.5);
319     FileOffset->SetMarkerStyle(21);
320     FileOffset->Draw("9p");*/
321 pam-rm2 1.1
322     //-------to save---------------------------///
323     oss.str("");
324     oss1.str("");
325    
326     oss << outDir.Data() << filename.Data();
327     oss1 << outDir.Data() << filename.Data();
328    
329     oss << "_HeaderScan_1." << format.Data();
330     oss1 << "_HeaderScan_2." << format.Data();
331    
332     finalCanv1->SaveAs(oss.str().c_str());
333 pam-rm2 1.1.1.2 //finalCanv2->SaveAs(oss1.str().c_str());
334 pam-rm2 1.1
335     file->Close();
336    
337     }
338    
339    
340    
341     int main(int argc, char* argv[]){
342     TString path;
343     TString outDir = "./";
344     TString format = "jpg";
345     if (argc < 2){
346     printf("You have to insert at least the file to analyze \n");
347     printf("Try '--help' for more information. \n");
348     exit(1);
349     }
350     if (!strcmp(argv[1], "--help")){
351     printf( "Usage: HeaderScan FILE [OPTION] \n");
352     printf( "\t --help Print this help and exit \n");
353     printf( "\t -outDir[path] Path where to put the output [default ./] \n");
354     printf( "\t -format[jpg|ps|gif] Format for output files [default 'jpg'] \n");
355     exit(1);
356     }
357     path=argv[1];
358 pam-rm2 1.1.1.2 for (int i = 2; i < argc; i++){
359 pam-rm2 1.1 if (!strcmp(argv[i], "-outDir")){
360     if (++i >= argc){
361     printf( "-outDir needs arguments. \n");
362     printf( "Try '--help' for more information. \n");
363     exit(1);
364     }
365     else{
366     outDir = argv[i];
367     continue;
368     }
369 pam-rm2 1.1.1.2 }
370 pam-rm2 1.1 if (!strcmp(argv[i], "-format")){
371     if (++i >= argc){
372     printf( "-format needs arguments. \n");
373     printf( "Try '--help' for more information. \n");
374     exit(1);
375     }
376     else{
377     format = argv[i];
378     continue;
379     }
380     }
381     }
382     HeaderScan(argv[1], outDir, format);
383     }

  ViewVC Help
Powered by ViewVC 1.1.23