]>
Commit | Line | Data |
---|---|---|
1 | /* | |
2 | TPC DA for online calibration | |
3 | ||
4 | Contact: Haavard.Helstrup@cern.ch | |
5 | Link: | |
6 | Run Type: PEDESTAL | |
7 | DA Type: LDC | |
8 | Number of events needed: 100 | |
9 | Input Files: | |
10 | Output Files: tpcPedestal.root, to be exported to the DAQ FXS | |
11 | fileId: pedestals | |
12 | Trigger types used: CALIBRATION_EVENT | |
13 | ||
14 | */ | |
15 | ||
16 | /* | |
17 | ||
18 | TPCda_pedestal.cxx - calibration algorithm for TPC pedestal runs | |
19 | ||
20 | 10/06/2007 sylvain.chapeland@cern.ch : first version - clean skeleton based on DAQ DA case1 | |
21 | 19/10/2007 christian.lippmann@cern.ch : Possibility to write output to ASCII file | |
22 | 24/10/2007 christian.lippmann@cern.ch : Including pedestal calibration for time bins | |
23 | 23/11/2007 christian.lippmann@cern.ch : Fix in order to avoid streamer problems in case of | |
24 | invalid ROOTSTYS. The famous magic line provided by Rene. | |
25 | 28/11/2007 christian.lippmann@cern.ch : TPC mapping file is read from DaqDetDB | |
26 | ||
27 | contact: marian.ivanov@cern.ch | |
28 | ||
29 | This process reads RAW data from the files provided as command line arguments | |
30 | and save results in a file (named from RESULT_FILE define - see below). | |
31 | ||
32 | */ | |
33 | ||
34 | #define RESULT_FILE "tpcPedestal.root" | |
35 | #define FILE_ID "pedestals" | |
36 | #define MAPPING_FILE "tpcMapping.root" | |
37 | #define AliDebugLevel() -1 | |
38 | ||
39 | extern "C" { | |
40 | #include <daqDA.h> | |
41 | } | |
42 | #include "event.h" | |
43 | #include "monitor.h" | |
44 | ||
45 | #include "stdio.h" | |
46 | #include "stdlib.h" | |
47 | #include <fstream> | |
48 | ||
49 | // | |
50 | //Root includes | |
51 | // | |
52 | #include "TFile.h" | |
53 | #include "TArrayF.h" | |
54 | #include "TROOT.h" | |
55 | #include "TPluginManager.h" | |
56 | ||
57 | // | |
58 | //AliRoot includes | |
59 | // | |
60 | #include "AliRawReader.h" | |
61 | #include "AliRawReaderDate.h" | |
62 | #include "AliTPCmapper.h" | |
63 | #include "AliTPCRawStream.h" | |
64 | #include "AliTPCROC.h" | |
65 | #include "AliTPCCalROC.h" | |
66 | #include "AliTPCCalPad.h" | |
67 | #include "AliMathBase.h" | |
68 | #include "TTreeStream.h" | |
69 | #include "AliLog.h" | |
70 | #include "TSystem.h" | |
71 | ||
72 | // | |
73 | // TPC calibration algorithm includes | |
74 | // | |
75 | #include "AliTPCCalibPedestal.h" | |
76 | ||
77 | /* | |
78 | Main routine, TPC pedestal detector algorithm to be run on TPC LDC | |
79 | Arguments: list of DATE raw data files | |
80 | */ | |
81 | ||
82 | int main(int argc, char **argv) { | |
83 | // | |
84 | // Main for TPC pedestal detector algorithm | |
85 | // | |
86 | ||
87 | AliLog::SetClassDebugLevel("AliTPCRawStream",-5); | |
88 | AliLog::SetClassDebugLevel("AliRawReaderDate",-5); | |
89 | AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5); | |
90 | AliLog::SetModuleDebugLevel("RAW",-5); | |
91 | ||
92 | Bool_t timeAnalysis = kTRUE; | |
93 | ||
94 | if (argc<2) { | |
95 | printf("Wrong number of arguments\n"); | |
96 | return -1; | |
97 | } | |
98 | ||
99 | /* magic line */ | |
100 | gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo", | |
101 | "*", | |
102 | "TStreamerInfo", | |
103 | "RIO", | |
104 | "TStreamerInfo()"); | |
105 | int i, status; | |
106 | ||
107 | /* log start of process */ | |
108 | printf("TPC DA started - %s\n",__FILE__); | |
109 | ||
110 | /* declare monitoring program */ | |
111 | status=monitorDeclareMp( __FILE__ ); | |
112 | if (status!=0) { | |
113 | printf("monitorDeclareMp() failed : %s\n",monitorDecodeError(status)); | |
114 | return -1; | |
115 | } | |
116 | ||
117 | AliTPCmapper *mapping = 0; // The TPC mapping | |
118 | ||
119 | if (!mapping){ | |
120 | /* copy locally the mapping file from daq detector config db */ | |
121 | status = daqDA_DB_getFile(MAPPING_FILE,"./tpcMapping.root"); | |
122 | if (status) { | |
123 | printf("Failed to get mapping file (%s) from DAQdetDB, status=%d\n", MAPPING_FILE, status); | |
124 | printf("Continue anyway ... maybe it works?\n"); // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!! | |
125 | return -1; // temporarily uncommented for testing on pcald47 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!! | |
126 | } | |
127 | ||
128 | /* open the mapping file and retrieve mapping object */ | |
129 | TFile *fileMapping = new TFile(MAPPING_FILE, "read"); | |
130 | mapping = (AliTPCmapper*) fileMapping->Get("tpcMapping"); | |
131 | delete fileMapping; | |
132 | } | |
133 | ||
134 | if (mapping == 0) { | |
135 | printf("Failed to get mapping object from %s. ...\n", MAPPING_FILE); | |
136 | //return -1; | |
137 | } else { | |
138 | printf("Got mapping object from %s\n", MAPPING_FILE); | |
139 | } | |
140 | ||
141 | AliTPCCalibPedestal calibPedestal; // pedestal and noise calibration | |
142 | calibPedestal.SetRangeTime(60,940); // set time bin range | |
143 | calibPedestal.SetTimeAnalysis(timeAnalysis); // pedestal(t) calibration | |
144 | if (mapping){ | |
145 | calibPedestal.SetAltroMapping(mapping->GetAltroMapping()); // Use altro mapping we got from daqDetDb | |
146 | } | |
147 | /* loop over RAW data files */ | |
148 | int nevents=0; | |
149 | for ( i=1; i<argc; i++ ) { | |
150 | ||
151 | /* define data source : this is argument i */ | |
152 | printf("Processing file %s\n", argv[i]); | |
153 | status=monitorSetDataSource( argv[i] ); | |
154 | if (status!=0) { | |
155 | printf("monitorSetDataSource() failed. Error=%s. Exiting ...\n", monitorDecodeError(status)); | |
156 | return -1; | |
157 | } | |
158 | ||
159 | /* read until EOF */ | |
160 | for ( ; ; ) { | |
161 | struct eventHeaderStruct *event; | |
162 | ||
163 | /* check shutdown condition */ | |
164 | if (daqDA_checkShutdown()) {break;} | |
165 | ||
166 | /* get next event (blocking call until timeout) */ | |
167 | status=monitorGetEventDynamic((void **)&event); | |
168 | if (status==MON_ERR_EOF) { | |
169 | printf ("End of File %d (%s) detected\n", i, argv[i]); | |
170 | break; /* end of monitoring file has been reached */ | |
171 | } | |
172 | if (status!=0) { | |
173 | printf("monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status)); | |
174 | break; | |
175 | } | |
176 | ||
177 | /* skip start/end of run events */ | |
178 | if ( (event->eventType != physicsEvent) && (event->eventType != calibrationEvent) ) | |
179 | continue; | |
180 | ||
181 | /* retry if got no event */ | |
182 | if (event==NULL) | |
183 | continue; | |
184 | ||
185 | nevents++; | |
186 | ||
187 | // Pedestal calibration | |
188 | AliRawReader *rawReader = new AliRawReaderDate((void*)event); | |
189 | calibPedestal.ProcessEvent(rawReader); | |
190 | //calibPedestal.ProcessEventFast(rawReader); // fast data reader | |
191 | delete rawReader; | |
192 | ||
193 | /* free resources */ | |
194 | free(event); | |
195 | } | |
196 | } | |
197 | ||
198 | // | |
199 | // Analyse pedestals and write them to rootfile | |
200 | // | |
201 | ||
202 | calibPedestal.Analyse(); | |
203 | calibPedestal.AnalyseTime(nevents); | |
204 | printf ("%d physics/calibration events processed.\n",nevents); | |
205 | ||
206 | TFile *fileTPC = new TFile(RESULT_FILE, "recreate"); | |
207 | calibPedestal.Write("tpcCalibPedestal"); | |
208 | delete fileTPC; | |
209 | printf("Wrote %s.\n",RESULT_FILE); | |
210 | ||
211 | /* store the result file on FES */ | |
212 | ||
213 | status=daqDA_FES_storeFile(RESULT_FILE,FILE_ID); | |
214 | if (status) { | |
215 | status = -2; | |
216 | } | |
217 | ||
218 | ||
219 | // | |
220 | // Now prepare ASCII files for local ALTRO configuration through DDL. | |
221 | // | |
222 | ||
223 | ofstream pedfile; | |
224 | ofstream noisefile; | |
225 | ofstream pedmemfile; | |
226 | char filename[255]; | |
227 | sprintf(filename,"tpcPedestals.data"); | |
228 | pedfile.open(filename); | |
229 | sprintf(filename,"tpcNoise.data"); | |
230 | noisefile.open(filename); | |
231 | sprintf(filename,"tpcPedestalMem.data"); | |
232 | pedmemfile.open(filename); | |
233 | ||
234 | TArrayF **timePed = calibPedestal.GetTimePedestals(); // pedestal values for each time bin | |
235 | ||
236 | Int_t ctr_channel = 0; | |
237 | Int_t ctr_altro = 0; | |
238 | Int_t ctr_pattern = 0; | |
239 | ||
240 | pedfile << 10 << std::endl; // mark file to contain PEDESTALS per channel | |
241 | noisefile << 11 << std::endl; // mark file to contain NOISE per altro | |
242 | pedmemfile << 12 << std::endl; // mark file to contain PEDESTALs per time bin | |
243 | ||
244 | for ( Int_t roc = 0; roc < 72; roc++ ) { | |
245 | if ( !calibPedestal.GetCalRocPedestal(roc) ) continue; | |
246 | Int_t side = mapping->GetSideFromRoc(roc); | |
247 | Int_t sector = mapping->GetSectorFromRoc(roc); | |
248 | //printf("Analysing ROC %d (side %d, sector %d) ...\n", roc, side, sector); | |
249 | Int_t nru = mapping->IsIROC(roc) ? 2 : 4; | |
250 | for ( int rcu = 0; rcu < nru; rcu++ ) { | |
251 | Int_t patch = mapping->IsIROC(roc) ? rcu : rcu+2; | |
252 | for ( int branch = 0; branch < 2; branch++ ) { | |
253 | for ( int fec = 0; fec < mapping->GetNfec(patch, branch); fec++ ) { | |
254 | for ( int altro = 0; altro < 8; altro++ ) { | |
255 | Float_t rms = 0.; | |
256 | Float_t ctr = 0.; | |
257 | for ( int channel = 0; channel < 16; channel++ ) { | |
258 | Int_t hwadd = mapping->CodeHWAddress(branch, fec, altro, channel); | |
259 | Int_t row = mapping->GetPadRow(patch, hwadd); // row in a ROC | |
260 | Int_t globalrow = mapping->GetGlobalPadRow(patch, hwadd); // row in full sector | |
261 | Int_t pad = mapping->GetPad(patch, hwadd); | |
262 | Float_t ped = calibPedestal.GetCalRocPedestal(roc)->GetValue(row,pad); | |
263 | // fixed pedestal | |
264 | if ( ped > 1.e-10 ) { | |
265 | pedfile << ctr_channel << "\t" << side << "\t" << sector << "\t" << patch << "\t" | |
266 | << hwadd << "\t" << ped << std::endl; | |
267 | ctr_channel++; | |
268 | } | |
269 | // pedestal(t) | |
270 | if ( timePed && fabs(timePed[globalrow][pad].GetSum()) > 1e-10 ) { | |
271 | pedmemfile << ctr_pattern << "\t" << side << "\t" << sector << "\t" << patch | |
272 | << "\t" << hwadd; | |
273 | for ( Int_t timebin = 0; timebin < 1024; timebin++ ) | |
274 | pedmemfile << "\t" << timePed[globalrow][pad].At(timebin); | |
275 | pedmemfile << std::endl; | |
276 | ctr_pattern++; | |
277 | } | |
278 | // rms=noise | |
279 | Float_t rms2 = calibPedestal.GetCalRocRMS(roc)->GetValue(row,pad); | |
280 | if ( rms2 > 1.e-10 ) { rms += rms2; ctr += 1.; } | |
281 | } // end channel for loop | |
282 | // noise data (rms) averaged over all channels in this ALTRO. | |
283 | Int_t hwadd = mapping->CodeHWAddress(branch, fec, altro, 0); | |
284 | if ( ctr > 1.e-10 ) { | |
285 | noisefile << ctr_altro << "\t" << side << "\t" << sector << "\t" << patch << "\t" | |
286 | << hwadd << "\t" << rms/ctr << std::endl; | |
287 | ctr_altro++; | |
288 | } | |
289 | } // end altro for loop | |
290 | } // end fec for loop | |
291 | } // end branch for loop | |
292 | } // end rcu for loop | |
293 | } // end roc loop | |
294 | ||
295 | pedfile.close(); | |
296 | noisefile.close(); | |
297 | pedmemfile.close(); | |
298 | printf("Wrote ASCII files.\n"); | |
299 | ||
300 | ||
301 | return status; | |
302 | } |