Call SetActiveSectors after the next event is loaded
[u/mrichter/AliRoot.git] / T0 / T0Laserda.cxx
CommitLineData
f16a53cf 1/*
2TOF DA for online calibration
3
4Contact: Michal.Oledzki@cern.ch
5Link: www.if.pw.edu.pl:/~oledzki
6Run Type: PHYSICS
7DA Type: MON
8Number of events needed: depending on the run, being run-level
9Input Files:
10Output Files: daLaser.root, to be exported to the DAQ FXS
11Trigger types used: PHYSICS_EVENT
12
13*/
99194e07 14
f16a53cf 15#define FILE_OUT "daLaser.root"
99194e07 16
f16a53cf 17#include <daqDA.h>
99194e07 18#include <event.h>
19#include <monitor.h>
20
21#include <Riostream.h>
22#include <stdio.h>
23#include <stdlib.h>
24
25//AliRoot
26#include <AliRawReaderDate.h>
27#include <AliRawReader.h>
28#include <AliT0RawReader.h>
99194e07 29
30//ROOT
0acfe114 31#include "TROOT.h"
32#include "TPluginManager.h"
33#include "TFile.h"
99194e07 34#include "TFile.h"
35#include "TKey.h"
36#include "TH2S.h"
37#include "TObject.h"
38#include "TBenchmark.h"
39#include "TRandom.h"
40#include "TMath.h"
41#include "TCanvas.h"
42#include "TString.h"
43#include "TH1.h"
44#include "TF1.h"
45#include "TSpectrum.h"
46#include "TVirtualFitter.h"
47#include "TProfile.h"
48
49/* Main routine
50 Arguments:
51 1- monitoring data source
52*/
53int main(int argc, char **argv) {
54 int status;
55
0acfe114 56 /* magic line */
57 gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
58 "*",
59 "TStreamerInfo",
60 "RIO",
61 "TStreamerInfo()");
99194e07 62
63 if (argc!=2) {
64 printf("Wrong number of arguments\n");
65 return -1;
66 }
67
68
69 /* define data source : this is argument 1 */
70 status=monitorSetDataSource( argv[1] );
71 if (status!=0) {
72 printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
73 return -1;
74 }
75
76
77 /* declare monitoring program */
78 status=monitorDeclareMp( __FILE__ );
79 if (status!=0) {
80 printf("monitorDeclareMp() failed : %s\n",monitorDecodeError(status));
81 return -1;
82 }
83
84
85 /* define wait event timeout - 1s max */
86 monitorSetNowait();
87 monitorSetNoWaitNetworkTimeout(1000);
88
89
90 /* log start of process */
91 printf("T0 monitoring program started\n");
92
93 // Allocation of histograms - start
94 TH1F *hCFD[24];
0acfe114 95 TH2F *hCFD_QTC[24];
96 TH2F *hCFD_LED[24];
99194e07 97
98 for(Int_t ic=0; ic<24; ic++) {
0acfe114 99 hCFD_QTC[ic] = new TH2F(Form("CFD_QTC%d",ic+1),"CFD_QTC",800,0.5,8000.5,200,10000.5,20000.5);
100 hCFD_LED[ic] = new TH2F(Form("CFD_LED%d",ic+1),"CFD_LED",100,-500.0,500.0,100,14600.0,14700.0);
99194e07 101 if(ic<12){
0acfe114 102 hCFD[ic] = new TH1F(Form("T0_C_%d_CFD",ic+1),"CFD", 3000,0., 30000.);
99194e07 103 }
104 else{
0acfe114 105 hCFD[ic] = new TH1F(Form("T0_A_%d_CFD",ic-11),"CFD", 3000,0., 30000.);
99194e07 106 }
107 }
108
109 // Allocation of histograms - end
110
111 Int_t iev=0;
112 /* main loop (infinite) */
113 for(;;) {
114 struct eventHeaderStruct *event;
115 eventTypeType eventT;
116
117 /* check shutdown condition */
118 if (daqDA_checkShutdown()) {break;}
119
120 /* get next event (blocking call until timeout) */
121 status=monitorGetEventDynamic((void **)&event);
122 if (status==(int)MON_ERR_EOF) {
123 printf ("End of File detected\n");
124 break; /* end of monitoring file has been reached */
125 }
126
127 if (status!=0) {
128 printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
129 break;
130 }
131
132 /* retry if got no event */
133 if (event==NULL) {
134 continue;
135 }
136
137 iev++;
138
139 /* use event - here, just write event id to result file */
140 eventT=event->eventType;
141
142 switch (event->eventType){
143
144 case START_OF_RUN:
145 break;
146
147 case END_OF_RUN:
148 break;
149
150// case CALIBRATION_EVENT:
151 case PHYSICS_EVENT:
152 printf(" event number = %i \n",iev);
153
154
155 // Initalize raw-data reading and decoding
156 AliRawReader *reader = new AliRawReaderDate((void*)event);
157
158 // Enable the following two lines in case of real-data
159 // reader->LoadEquipmentIdsMap("T0map.txt");
160 // reader->RequireHeader(kFALSE);
161 // reader->RequireHeader(kTRUE);
162 AliT0RawReader *start = new AliT0RawReader(reader, kTRUE);
163
164 // Read raw data
165 Int_t allData[105][5];
166 for(Int_t i0=0;i0<105;i0++)
167 for(Int_t j0=0;j0<5;j0++)
168 allData[i0][j0] = 0;
169
170 if(start->Next())
171 for (Int_t i=0; i<105; i++) {
172 for(Int_t iHit=0;iHit<5;iHit++){
173 allData[i][iHit]= start->GetData(i,iHit);
174 }
175 }
176 else
177 printf("No T0 data found!!\n");
178
179 // Fill the histograms
180
181 for (Int_t ik = 0; ik<24; ik+=2)
182 for (Int_t iHt=0; iHt<5; iHt++){
183 Int_t cc = ik/2;
184 if((allData[cc+1][iHt]-allData[0][0])>0){
185 hCFD[cc]->Fill(allData[cc+1][iHt]-allData[0][0]);
186 }
187 if(allData[ik+25][iHt]!=0 && allData[ik+26][iHt]!=0 && allData[cc+1][iHt]!=0){
188 hCFD_QTC[cc]->Fill((allData[ik+25][iHt]-allData[ik+26][iHt]) , (allData[cc+1][iHt]-allData[0][0]));
189 }
190 if(allData[cc+13][iHt]!=0 && allData[cc+1][iHt]!=0){
191 hCFD_LED[cc]->Fill(allData[cc+13][iHt]-allData[cc+1][iHt],allData[cc+1][iHt]-allData[0][0]);
192 }
193 }
194
195 for (Int_t ik = 24; ik<48; ik+=2)
196 for (Int_t iHt=0; iHt<5; iHt++){
197 Int_t cc = ik/2;
198 if((allData[cc+45][iHt]-allData[0][0])>0){
199 hCFD[cc]->Fill(allData[cc+45][iHt]-allData[0][0]);
200 }
201
202 if(allData[ik+57][iHt]!=0 && allData[ik+58][iHt]!=0){
203 if(allData[cc+45][iHt]!=0)
204 hCFD_QTC[cc]->Fill(allData[ik+57][iHt]-allData[ik+58][iHt],allData[cc+45][iHt]-allData[0][0]);
205 }
206 if(allData[cc+57][iHt]!=0 && allData[cc+45][iHt]!=0){
207 hCFD_LED[cc]->Fill(allData[cc+57][iHt]-allData[cc+45][iHt],allData[cc+45][iHt]-allData[0][0]);
208 }
209 }
210
211
212 delete start;
213 start = 0x0;
214 reader->Reset();
215 // End of fill histograms
216
217 }
218
219
220
221 /* free resources */
222 free(event);
223
224 /* exit when last event received, no need to wait for TERM signal */
225 if (eventT==END_OF_RUN) {
226 printf("EOR event detected\n");
227 break;
228 }
229 }
230printf("After loop, before writing histos\n");
231 // write a file with the histograms
f16a53cf 232 TFile *hist = new TFile(FILE_OUT,"RECREATE");
99194e07 233
234 for(Int_t j=0;j<24;j++){
235 hCFD_QTC[j]->Write();
236 hCFD_LED[j]->Write();
237 hCFD[j]->Write();
238 }
239 hist->Close();
240 delete hist;
241
242 return status;
243}
244