]>
Commit | Line | Data |
---|---|---|
1 | /* | |
2 | ||
3 | TOF DA for online calibration | |
4 | ||
5 | Contact: Chiara.Zampolli@bo.infn.it | |
6 | Link: www.bo.infn.it/~zampolli | |
7 | Run Type: PHYSICS | |
8 | DA Type: MON | |
9 | Number of events needed: depending on the run, being run-level | |
10 | Input Files: TOFdaTotal.root, to be updated if existing | |
11 | Output Files: TOFdaRun.root, TOFdaTotal.root, both to be exported to the DAQ FXS | |
12 | Trigger types used: PHYSICS_EVENT | |
13 | ||
14 | */ | |
15 | ||
16 | #define FILE_TOTAL "TOFdaTotal.root" | |
17 | #define FILE_RUN "TOFdaRun.root" | |
18 | ||
19 | // DATE | |
20 | #include <daqDA.h> | |
21 | #include <event.h> | |
22 | #include <monitor.h> | |
23 | ||
24 | //AliRoot | |
25 | #include <AliTOFRawStream.h> | |
26 | #include <AliRawReaderDate.h> | |
27 | #include <AliRawReader.h> | |
28 | #include <AliTOFGeometry.h> | |
29 | #include <AliT0RawReader.h> | |
30 | #include <AliDAQ.h> | |
31 | #include <AliTOFHitData.h> | |
32 | #include <AliTOFHitDataBuffer.h> | |
33 | ||
34 | //ROOT | |
35 | #include <TFile.h> | |
36 | #include <TKey.h> | |
37 | #include <TH2S.h> | |
38 | #include <TObject.h> | |
39 | #include <TMath.h> | |
40 | #include <TSystem.h> | |
41 | ||
42 | ||
43 | /* Main routine | |
44 | Arguments: | |
45 | 1- monitoring data source | |
46 | */ | |
47 | ||
48 | int main(int argc, char **argv) { | |
49 | ||
50 | AliTOFGeometry * geom = new AliTOFGeometry(); | |
51 | ||
52 | static const Int_t size = AliTOFGeometry::NPadXSector()*AliTOFGeometry::NSectors(); | |
53 | static const Int_t nbins = 500; | |
54 | static const Int_t binmin = -20; | |
55 | const Float_t c = 2.99792458E10; //speed of light | |
56 | TH1F::AddDirectory(0); | |
57 | TH2S * htofPartial = new TH2S("htof","histo with delays", size,-0.5,size*1.-0.5,nbins,binmin-0.5,nbins*1.+binmin-0.5); | |
58 | ||
59 | // decoding the events | |
60 | ||
61 | int status; | |
62 | if (argc!=2) { | |
63 | printf("Wrong number of arguments\n"); | |
64 | return -1; | |
65 | } | |
66 | ||
67 | /* define data source : this is argument 1 */ | |
68 | status=monitorSetDataSource( argv[1] ); | |
69 | if (status!=0) { | |
70 | printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status)); | |
71 | return -1; | |
72 | } | |
73 | ||
74 | /* declare monitoring program */ | |
75 | status=monitorDeclareMp( __FILE__ ); | |
76 | if (status!=0) { | |
77 | printf("monitorDeclareMp() failed : %s\n",monitorDecodeError(status)); | |
78 | return -1; | |
79 | } | |
80 | ||
81 | /* define wait event timeout - 1s max */ | |
82 | monitorSetNowait(); | |
83 | monitorSetNoWaitNetworkTimeout(1000); | |
84 | ||
85 | /* log start of process */ | |
86 | printf("TOF DA started\n"); | |
87 | ||
88 | /* init some counters */ | |
89 | int nevents_physics=0; | |
90 | int nevents_total=0; | |
91 | ||
92 | struct eventHeaderStruct *event; | |
93 | eventTypeType eventT; | |
94 | Int_t iev=0; | |
95 | ||
96 | /* main loop (infinite) */ | |
97 | for(;;) { | |
98 | ||
99 | /* check shutdown condition */ | |
100 | if (daqDA_checkShutdown()) {break;} | |
101 | ||
102 | /* get next event (blocking call until timeout) */ | |
103 | status=monitorGetEventDynamic((void **)&event); | |
104 | if (status==MON_ERR_EOF) { | |
105 | printf ("End of File detected\n"); | |
106 | break; /* end of monitoring file has been reached */ | |
107 | } | |
108 | ||
109 | if (status!=0) { | |
110 | printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status)); | |
111 | break; | |
112 | } | |
113 | ||
114 | /* retry if got no event */ | |
115 | if (event==NULL) { | |
116 | continue; | |
117 | } | |
118 | ||
119 | iev++; | |
120 | ||
121 | /* use event - here, just write event id to result file */ | |
122 | nevents_total++; | |
123 | eventT=event->eventType; | |
124 | switch (event->eventType){ | |
125 | ||
126 | /* START OF RUN */ | |
127 | case START_OF_RUN: | |
128 | break; | |
129 | /* END START OF RUN */ | |
130 | ||
131 | /* END OF RUN */ | |
132 | case END_OF_RUN: | |
133 | break; | |
134 | ||
135 | case PHYSICS_EVENT: | |
136 | nevents_physics++; | |
137 | AliRawReader *rawReader = new AliRawReaderDate((void*)event); | |
138 | //rawReader->RequireHeader(kFALSE); | |
139 | ||
140 | //T0 event | |
141 | Int_t meantime = 0; | |
142 | AliT0RawReader *rawReaderT0 = new AliT0RawReader(rawReader,kTRUE); | |
143 | if (!rawReaderT0->Next()) { | |
144 | printf("T0: no raw data found!\n"); | |
145 | } | |
146 | else { | |
147 | /* | |
148 | Int_t allData[105][5]; | |
149 | for (Int_t i=0; i<105; i++) { | |
150 | allData[i][0]=rawReaderT0->GetData(i,0); | |
151 | } | |
152 | meantime = allData[49][0]; | |
153 | */ | |
154 | meantime = rawReaderT0->GetData(49,0); | |
155 | // printf("time zero (ns) = %i (%f) \n", meantime, (meantime*24.4-200)*1E-3); // debugging purpose | |
156 | } | |
157 | ||
158 | delete rawReaderT0; | |
159 | rawReaderT0 = 0x0; | |
160 | rawReader->Reset(); | |
161 | ||
162 | //TOF event | |
163 | Int_t dummy = -1; | |
164 | Int_t Volume[5]; | |
165 | AliTOFHitData *HitData; | |
166 | AliTOFHitDataBuffer *DataBuffer; | |
167 | AliTOFHitDataBuffer *PackedDataBuffer; | |
168 | AliTOFRawStream *rawStreamTOF = new AliTOFRawStream(rawReader); | |
169 | // rawReader->ReadHeader(); | |
170 | rawStreamTOF->ResetBuffers(); | |
171 | rawStreamTOF->DecodeDDL(0, AliDAQ::NumberOfDdls("TOF") - 1,0); | |
172 | Int_t nPDBEntriesToT = 0; | |
173 | Int_t nDBEntriesToT = 0; | |
174 | for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls("TOF"); iDDL++){ | |
175 | ||
176 | /* read decoded data */ | |
177 | DataBuffer = rawStreamTOF->GetDataBuffer(iDDL); | |
178 | PackedDataBuffer = rawStreamTOF->GetPackedDataBuffer(iDDL); | |
179 | ||
180 | /* get buffer entries */ | |
181 | Int_t nDBEntries = DataBuffer->GetEntries(); | |
182 | Int_t nPDBEntries = PackedDataBuffer->GetEntries(); | |
183 | nPDBEntriesToT+=nPDBEntries; | |
184 | nDBEntriesToT+=nDBEntries; | |
185 | ||
186 | //for (Int_t iHit = 0; iHit < nDBEntries; iHit++){ | |
187 | // HitData = DataBuffer->GetHit(iHit); | |
188 | /* store volume information */ | |
189 | // rawStreamTOF->EquipmentId2VolumeId(HitData, Volume); | |
190 | //} | |
191 | /* reset buffer */ | |
192 | DataBuffer->Reset(); | |
193 | ||
194 | /* read data buffer hits */ | |
195 | for (Int_t iHit = 0; iHit < nPDBEntries; iHit++){ | |
196 | HitData = PackedDataBuffer->GetHit(iHit); | |
197 | /* add volume information */ | |
198 | HitData->SetDDLID(iDDL); | |
199 | rawStreamTOF->EquipmentId2VolumeId(HitData, Volume); | |
200 | if (Volume[0]==-1 || | |
201 | Volume[1]==-1 || | |
202 | Volume[2]==-1 || | |
203 | Volume[3]==-1 || | |
204 | Volume[4]==-1) continue; | |
205 | else { | |
206 | dummy = Volume[3]; | |
207 | Volume[3] = Volume[4]; | |
208 | Volume[4] = dummy; | |
209 | Int_t tof = (Int_t)((Double_t)HitData->GetTime()*1E3/AliTOFGeometry::TdcBinWidth()); | |
210 | Int_t index = geom->GetIndex(Volume); | |
211 | Float_t pos[3]; | |
212 | geom->GetPosPar(Volume,pos); | |
213 | Float_t texp=TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]+pos[2]*pos[2])/c*1E9; //expected time in ns | |
214 | Float_t texpBin=(texp*1E3-32)/AliTOFGeometry::TdcBinWidth(); //expected time in number of TDC bin | |
215 | Int_t deltabin = tof-TMath::Nint(texpBin); //to be used with real data; rounding expected time to Int_t | |
216 | htofPartial->Fill(index,deltabin); //channel index start from 0, bin index from 1 | |
217 | //debugging printings | |
218 | //printf("sector %i, plate %i, strip %i, padz %i, padx %i \n",Volume[0],Volume[1],Volume[2],Volume[3],Volume[4]); | |
219 | //printf("pos x = %f, pos y = %f, pos z = %f \n",pos[0],pos[1],pos[2]); | |
220 | //printf ("expected time = %f (ns)\n",texp); | |
221 | //printf ("expected time bin = %f (TDC bin)\n",texpBin); | |
222 | //printf ("measured time bin = %i (TDC bin) with %f (ns) and ACQ bit = %i \n",tof, HitData->GetTime(), HitData->GetACQ()); | |
223 | // printf("index = %i, deltabin = %i , filling index = %i, and bin = % i\n",index, deltabin, index, deltabin); | |
224 | ||
225 | } | |
226 | /* reset buffer */ | |
227 | PackedDataBuffer->Reset(); | |
228 | } | |
229 | } | |
230 | //printf(" Packed Hit Buffer Entries = %i \n",nPDBEntriesToT); | |
231 | //printf(" Hit Buffer Entries = %i \n",nDBEntriesToT); | |
232 | ||
233 | delete rawStreamTOF; | |
234 | rawStreamTOF = 0x0; | |
235 | ||
236 | delete rawReader; | |
237 | rawReader = 0x0; | |
238 | } | |
239 | ||
240 | /* free resources */ | |
241 | free(event); | |
242 | ||
243 | /* exit when last event received, no need to wait for TERM signal */ | |
244 | if (eventT==END_OF_RUN) { | |
245 | printf("EOR event detected\n"); | |
246 | break; | |
247 | } | |
248 | } | |
249 | ||
250 | delete geom; | |
251 | geom = 0x0; | |
252 | ||
253 | //write the Run level file | |
254 | TFile * fileRun = new TFile (FILE_RUN,"RECREATE"); | |
255 | htofPartial->Write(); | |
256 | fileRun->Close(); | |
257 | ||
258 | //write the Total file | |
259 | TH2S *htoftot = 0x0; | |
260 | TFile * filetot = 0x0; | |
261 | Bool_t isThere=kFALSE; | |
262 | const char *dirname = "./"; | |
263 | TString filename = FILE_TOTAL; | |
264 | if((gSystem->FindFile(dirname,filename))!=NULL){ | |
265 | isThere=kTRUE; | |
266 | printf("%s found \n",FILE_TOTAL); | |
267 | } | |
268 | if(isThere){ | |
269 | ||
270 | TFile * filetot1 = new TFile (FILE_TOTAL,"READ"); | |
271 | //look for the file | |
272 | if (!filetot1->IsZombie()){ | |
273 | printf("updating file %s \n",FILE_TOTAL); | |
274 | TIter next(filetot1->GetListOfKeys()); | |
275 | TKey *key; | |
276 | //look for the histogram | |
277 | while ((key=(TKey*)next())){ | |
278 | const char * namekey = key->GetName(); | |
279 | if (strcmp(namekey,"htoftot")==0){ | |
280 | printf(" histo found \n"); | |
281 | htoftot = (TH2S*) filetot1->Get("htoftot"); | |
282 | htoftot->AddDirectory(0); | |
283 | htoftot->Add(htofPartial); | |
284 | break; | |
285 | } | |
286 | } | |
287 | } | |
288 | filetot1->Close(); | |
289 | delete filetot1; | |
290 | filetot1=0x0; | |
291 | } | |
292 | else { | |
293 | printf(" no %s file found \n",FILE_TOTAL); | |
294 | htoftot = new TH2S(*htofPartial); | |
295 | htoftot->SetName("htoftot"); | |
296 | htoftot->AddDirectory(0); | |
297 | } | |
298 | ||
299 | filetot = new TFile (FILE_TOTAL,"RECREATE"); | |
300 | filetot->cd(); | |
301 | htoftot->Write(); | |
302 | filetot->Close(); | |
303 | ||
304 | delete fileRun; | |
305 | delete filetot; | |
306 | delete htofPartial; | |
307 | delete htoftot; | |
308 | ||
309 | fileRun = 0x0; | |
310 | filetot = 0x0; | |
311 | htofPartial = 0x0; | |
312 | htoftot = 0x0; | |
313 | ||
314 | /* write report */ | |
315 | printf("Run #%s, received %d physics events out of %d\n",getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total); | |
316 | ||
317 | status=0; | |
318 | ||
319 | /* export file to FXS */ | |
320 | if (daqDA_FES_storeFile(FILE_RUN, FILE_RUN)) { | |
321 | status=-2; | |
322 | } | |
323 | if (daqDA_FES_storeFile(FILE_TOTAL, FILE_TOTAL)) { | |
324 | status=-2; | |
325 | } | |
326 | ||
327 | ||
328 | return status; | |
329 | } | |
330 |