]>
Commit | Line | Data |
---|---|---|
29c170aa | 1 | /* |
2 | ||
3 | TOF DA for online calibration from pulser data | |
4 | ||
5 | Contact: Chiara.Zampolli@bo.infn.it | |
6 | Link: www.bo.infn.it/~zampolli | |
7fffa85b | 7 | Run Type: NOISE |
29c170aa | 8 | DA Type: LDC |
9 | Number of events needed: 10000 | |
10 | Input Files: TOF<nrun>.raw, where <nrun> is the run number | |
11 | Output Files: TOFoutPulserLDC_<LDCid>.root, where <LDCid> is the id of the LDC which is read (2 characters field, e.g. TOFoutPulserLDC_03.root),to be exported to the DAQ FXS | |
12 | Trigger types used: PHYSICS_EVENT (for the time being) | |
13 | ||
14 | */ | |
15 | ||
16 | // DATE | |
17 | #include "event.h" | |
18 | #include "monitor.h" | |
89b3ea8c | 19 | #include "daqDA.h" |
29c170aa | 20 | |
21 | #include <stdio.h> | |
22 | #include <stdlib.h> | |
23 | ||
24 | //AliRoot | |
25 | #include <AliTOFRawStream.h> | |
26 | #include <AliRawReaderDate.h> | |
27 | #include <AliRawReader.h> | |
28 | #include <AliTOFGeometry.h> | |
29 | #include <AliDAQ.h> | |
30 | #include <AliTOFHitData.h> | |
31 | #include <AliTOFHitDataBuffer.h> | |
32 | #include <AliTOFDecoder.h> | |
33 | ||
34 | //ROOT | |
35 | #include <TFile.h> | |
36 | #include <TKey.h> | |
37 | #include <TH1I.h> | |
38 | #include <TObject.h> | |
39 | #include <TMath.h> | |
40 | #include <TSystem.h> | |
3bc66139 | 41 | #include "TROOT.h" |
42 | #include "TPluginManager.h" | |
29c170aa | 43 | |
44 | /* Main routine | |
45 | Arguments: list of DATE raw data files | |
46 | */ | |
47 | int main(int argc, char **argv) { | |
48 | ||
3bc66139 | 49 | /* magic line from Rene */ |
50 | gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo", | |
51 | "*", | |
52 | "TStreamerInfo", | |
53 | "RIO", | |
54 | "TStreamerInfo()"); | |
55 | ||
29c170aa | 56 | AliTOFGeometry * geom = new AliTOFGeometry(); |
57 | ||
58 | static const Int_t size = AliTOFGeometry::NPadXSector()*AliTOFGeometry::NSectors(); | |
59 | TH1F::AddDirectory(0); | |
60 | TH1F * htofNoise = new TH1F("hTOFnoise","histo with signals on TOF during noise", size,-0.5,size-0.5); | |
61 | for (Int_t ibin =1;ibin<=size;ibin++){ | |
62 | htofNoise->SetBinContent(ibin,-1); | |
63 | } | |
64 | UInt_t ldcId=0; | |
65 | ||
66 | int status; | |
67 | ||
68 | /* log start of process */ | |
69 | printf("TOF DA started\n"); | |
70 | ||
71 | /* check that we got some arguments = list of files */ | |
72 | if (argc<2) { | |
73 | printf("Wrong number of arguments\n"); | |
74 | return -1; | |
75 | } | |
76 | ||
77 | /* open result file */ | |
78 | FILE *fp=NULL; | |
79 | fp=fopen("./result.txt","a"); | |
80 | if (fp==NULL) { | |
81 | printf("Failed to open file\n"); | |
82 | return -1; | |
83 | } | |
84 | ||
85 | /* init some counters */ | |
86 | int nevents_physics=0; | |
87 | int nevents_total=0; | |
88 | ||
89 | /* read the data files */ | |
90 | int n; | |
91 | for (n=1;n<argc;n++) { | |
92 | ||
93 | status=monitorSetDataSource( argv[n] ); | |
94 | if (status!=0) { | |
95 | printf("monitorSetDataSource() failed : %s\n",monitorDecodeError(status)); | |
96 | return -1; | |
97 | } | |
98 | ||
99 | /* read the file */ | |
100 | for(;;) { | |
101 | struct eventHeaderStruct *event; | |
102 | eventTypeType eventT; | |
103 | ||
104 | /* get next event */ | |
105 | status=monitorGetEventDynamic((void **)&event); | |
106 | if (status==MON_ERR_EOF) break; /* end of monitoring file has been reached */ | |
107 | if (status!=0) { | |
108 | printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status)); | |
109 | return -1; | |
110 | } | |
111 | ||
112 | /* retry if got no event */ | |
113 | if (event==NULL) { | |
114 | break; | |
115 | } | |
116 | ||
117 | /* use event - here, just write event id to result file */ | |
118 | eventT=event->eventType; | |
119 | ||
120 | if (eventT==PHYSICS_EVENT) { | |
121 | //printf ("event %i \n", nevents_physics); | |
122 | Int_t nPDBEntriesToT = 0; | |
123 | Int_t nDBEntriesToT = 0; | |
124 | AliTOFHitData *HitData; | |
125 | Int_t dummy = -1; | |
126 | Int_t Volume[5]; | |
127 | AliRawReader *rawReader = new AliRawReaderDate((void*)event); | |
128 | AliTOFRawStream *rawStreamTOF = new AliTOFRawStream(rawReader); | |
9512a964 | 129 | const AliRawDataHeader *currentCDH; |
29c170aa | 130 | AliTOFDecoder * decoderTOF = new AliTOFDecoder(); |
131 | AliTOFHitDataBuffer *DataBuffer[72]; | |
132 | AliTOFHitDataBuffer *PackedDataBuffer[72]; | |
133 | for (Int_t i=0;i<AliDAQ::NumberOfDdls("TOF");i++){ | |
134 | DataBuffer[i]=new AliTOFHitDataBuffer(); | |
135 | PackedDataBuffer[i]=new AliTOFHitDataBuffer(); | |
136 | } | |
137 | Int_t currentEquipment; | |
138 | Int_t currentDDL; | |
139 | UChar_t *data = 0x0; | |
140 | while(rawReader->ReadHeader()){ | |
141 | //printf ("ldcId = %i \n",ldcId); | |
142 | ldcId = rawReader->GetLDCId(); | |
143 | ||
144 | //memory leak prevention (actually data should be always 0x0 here) | |
145 | if (data != 0x0) delete [] data; | |
146 | ||
147 | //get equipment infos | |
148 | currentEquipment = rawReader->GetEquipmentId(); | |
149 | currentDDL = rawReader->GetDDLID(); | |
9512a964 | 150 | currentCDH = rawReader->GetDataHeader(); |
29c170aa | 151 | Int_t nchDDL = 0; |
152 | if (currentDDL%2==0) { | |
153 | nchDDL = 2160; | |
154 | } | |
155 | else { | |
156 | nchDDL = 2208; | |
157 | } | |
158 | Int_t * array = new Int_t[nchDDL]; | |
159 | decoderTOF->GetArrayDDL(array, currentDDL); | |
160 | ||
161 | for (Int_t i=0;i<nchDDL;i++){ | |
162 | if (htofNoise->GetBinContent(array[i]+1)<0) htofNoise->SetBinContent(array[i]+1,0); | |
163 | } | |
164 | ||
165 | //printf(" Equipment = %i, and DDL = %i \n", currentEquipment,currentDDL); | |
166 | const Int_t kDataSize = rawReader->GetDataSize(); | |
167 | const Int_t kDataWords = kDataSize / 4; | |
168 | data = new UChar_t[kDataSize]; | |
169 | decoderTOF->SetDataBuffer(DataBuffer[currentDDL]); | |
170 | decoderTOF->SetPackedDataBuffer(PackedDataBuffer[currentDDL]); | |
171 | //start decoding | |
172 | if (!rawReader->ReadNext(data, kDataSize)){ | |
173 | rawReader->AddMajorErrorLog(AliTOFRawStream::kDDLdataReading); | |
174 | printf("Error while reading DDL data. Go to next equipment \n"); | |
175 | delete [] data; | |
176 | data = 0x0; | |
177 | continue; | |
178 | } | |
9512a964 | 179 | if (decoderTOF->Decode((UInt_t *)data, kDataWords, currentCDH) == kTRUE ){ |
29c170aa | 180 | rawReader->AddMajorErrorLog(AliTOFRawStream::kDDLDecoder,Form("DDL # = %d",currentDDL)); |
181 | printf("Error while decoding DDL # %d: decoder returned with errors \n", currentDDL); | |
182 | } | |
183 | ||
184 | Int_t nDBEntries = DataBuffer[currentDDL]->GetEntries(); | |
185 | Int_t nPDBEntries = PackedDataBuffer[currentDDL]->GetEntries(); | |
186 | nPDBEntriesToT+=nPDBEntries; | |
187 | nDBEntriesToT+=nDBEntries; | |
188 | /* reset buffer */ | |
189 | DataBuffer[currentDDL]->Reset(); | |
190 | ||
191 | /* read data buffer hits */ | |
192 | for (Int_t iHit = 0; iHit < nPDBEntries; iHit++){ | |
193 | HitData = PackedDataBuffer[currentDDL]->GetHit(iHit); | |
194 | /* add volume information */ | |
195 | HitData->SetDDLID(currentDDL); | |
196 | rawStreamTOF->EquipmentId2VolumeId(HitData, Volume); | |
197 | if (Volume[0]==-1 || | |
198 | Volume[1]==-1 || | |
199 | Volume[2]==-1 || | |
200 | Volume[3]==-1 || | |
201 | Volume[4]==-1) continue; | |
202 | else { | |
203 | dummy = Volume[3]; | |
204 | Volume[3] = Volume[4]; | |
205 | Volume[4] = dummy; | |
206 | Int_t index = geom->GetIndex(Volume); | |
207 | Bool_t found =kFALSE; | |
208 | // to check array indexes | |
209 | /* | |
210 | for (Int_t j=0;j<nchDDL;j++){ | |
211 | if (index==array[j]) { | |
212 | found = kTRUE; | |
213 | break; | |
214 | } | |
215 | } | |
216 | printf ("index = %i, found = %i\n",index, (Int_t)found); | |
217 | */ | |
218 | //printf ("index = %i \n",index); | |
219 | htofNoise->Fill(index); //channel index start from 0, bin index from 1 | |
220 | //debugging printings | |
221 | //printf("sector %i, plate %i, strip %i, padz %i, padx %i \n",Volume[0],Volume[1],Volume[2],Volume[3],Volume[4]); | |
222 | } | |
223 | } | |
224 | /* reset buffer */ | |
225 | PackedDataBuffer[currentDDL]->Reset(); | |
226 | delete [] data; | |
227 | data = 0x0; | |
228 | delete [] array; | |
229 | } | |
230 | ||
231 | //printf(" Packed Hit Buffer Entries = %i \n",nPDBEntriesToT); | |
232 | //printf(" Hit Buffer Entries = %i \n",nDBEntriesToT); | |
233 | ||
234 | delete decoderTOF; | |
235 | decoderTOF=0x0; | |
236 | ||
237 | for (Int_t i=0;i<72;i++){ | |
238 | delete DataBuffer[i]; | |
239 | delete PackedDataBuffer[i]; | |
240 | } | |
241 | ||
242 | delete rawStreamTOF; | |
243 | rawStreamTOF = 0x0; | |
244 | ||
245 | delete rawReader; | |
246 | rawReader = 0x0; | |
247 | nevents_physics++; | |
248 | } | |
249 | nevents_total++; | |
250 | ||
251 | /* free resources */ | |
252 | free(event); | |
253 | } | |
254 | ||
255 | } | |
256 | delete geom; | |
257 | geom = 0x0; | |
258 | ||
259 | Float_t time = nevents_physics*200*1E-9; // acquisition time in s | |
260 | //printf(" Noise run lasted %f s \n",time); | |
261 | for (Int_t ibin =1;ibin<=size;ibin++){ | |
262 | Float_t cont = htofNoise->GetBinContent(ibin); | |
263 | if (cont!=-1) { | |
264 | //printf(" content = %f \n", cont); | |
265 | htofNoise->SetBinContent(ibin,cont/time); | |
266 | //printf(" scaled content = %f \n", cont/time); | |
267 | } | |
268 | } | |
269 | ||
270 | //write the Run level file | |
271 | char filename[100]; | |
272 | sprintf(filename,"TOFoutNoiseLDC_%02i.root",ldcId); | |
273 | TFile * fileRun = new TFile (filename,"RECREATE"); | |
274 | htofNoise->Write(); | |
275 | fileRun->Close(); | |
276 | ||
277 | /* write report */ | |
278 | fprintf(fp,"Run #%s, received %d physics events out of %d\n",getenv("DATE_RUN_NUMBER"),nevents_physics,nevents_total); | |
279 | ||
280 | /* close result file */ | |
281 | fclose(fp); | |
282 | ||
283 | status =0; | |
284 | ||
285 | /* store the result file on FES */ | |
98a837bf | 286 | status=daqDA_FES_storeFile(filename,"NOISE"); |
29c170aa | 287 | if (status) { |
288 | status = -2; | |
289 | } | |
290 | ||
291 | return status; | |
292 | } |