AliRoot Core  edcc906 (edcc906)
TPCCEda.cxx
Go to the documentation of this file.
1 /*
2 TPC DA for online calibration
3 
4 Contact: Haavard.Helstrup@cern.ch
5 Link:.
6 Run Type: PHYSICS STANDALONE DAQ
7 DA Type: MON
8 Number of events needed: 500
9 Input Files:.
10 Output Files: tpcCE.root, to be exported to the DAQ FXS
11 fileId: CE
12 Trigger types used: PHYSICS_EVENT
13 */
37 
38 #define RESULT_FILE "tpcCE.root"
39 #define FILE_ID "CE"
40 #define MAPPING_FILE "tpcMapping.root"
41 #define CONFIG_FILE "TPCCEda.conf"
42 
43 #include <daqDA.h>
44 #include "event.h"
45 #include "monitor.h"
46 #include <stdio.h>
47 #include <stdlib.h>
48 #include <vector>
49 //
50 //Root includes
51 //
52 #include <TFile.h>
53 #include "TROOT.h"
54 #include "TPluginManager.h"
55 #include "TSystem.h"
56 #include "TString.h"
57 #include "TObjString.h"
58 #include "TDatime.h"
59 #include "TStopwatch.h"
60 #include "TMap.h"
61 #include "TGraph.h"
62 #include "TMath.h"
63 //
64 //AliRoot includes
65 //
66 #include "AliRawReader.h"
67 #include "AliRawReaderDate.h"
68 #include "AliTPCmapper.h"
69 #include "AliTPCROC.h"
70 #include "AliTPCCalROC.h"
71 #include "AliTPCCalPad.h"
72 #include "AliMathBase.h"
73 #include "TTreeStream.h"
74 #include "AliLog.h"
75 #include "AliTPCConfigDA.h"
76 //
77 //AMORE
78 //
79 #include <AmoreDA.h>
80 //
81 // TPC calibration algorithm includes
82 //
83 #include "AliTPCCalibCE.h"
84 
85 
86 //functios, implementation below
87 void SendToAmoreDB(AliTPCCalibCE *calibCE, unsigned long32 runNb);
88 //for threaded processing
89 
90 
91 int main(int argc, char **argv) {
94 
95  // log start of process
96  printf("TPCCEda: DA started - %s\n",__FILE__);
97 
98  if (argc<2) {
99  printf("TPCCEda: Wrong number of arguments\n");
100  return -1;
101  }
102 
103  AliLog::SetClassDebugLevel("AliRawReaderDate",-5);
104  AliLog::SetClassDebugLevel("AliTPCAltroMapping",-5);
105  AliLog::SetModuleDebugLevel("RAW",-5);
106 
107  gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
108  "*",
109  "TStreamerInfo",
110  "RIO",
111  "TStreamerInfo()");
112 
113  /* declare monitoring program */
114  int i,status;
115  status=monitorDeclareMp( __FILE__ );
116  if (status!=0) {
117  printf("TPCCEda: monitorDeclareMp() failed : %s\n",monitorDecodeError(status));
118  return -1;
119  }
120  monitorSetNowait();
121  monitorSetNoWaitNetworkTimeout(1000);
122 
123  //variables
124  AliTPCmapper *mapping = 0; // The TPC mapping
125  char localfile[255];
126  unsigned long32 runNb=0; //run number
127 
128 
129  //
130  // DA configuration from configuration file
131  //
132  //retrieve configuration file
133  sprintf(localfile,"./%s",CONFIG_FILE);
134  status = daqDA_DB_getFile(CONFIG_FILE,localfile);
135  if (status) {
136  printf("TPCCEda: Failed to get configuration file (%s) from DAQdetDB, status=%d\n", CONFIG_FILE, status);
137  return -1;
138  }
139  AliTPCConfigDA config(CONFIG_FILE);
140  // check configuration options
141  TString laserTriggerName("C0LSR-ABCE-NOPF-CENT");
142  TString monitoringType("YES");
143  Int_t forceTriggerId=-1;
144  Int_t saveOption=2; // how to store the object. See AliTPCCalibCE::DumpToFile
145  Bool_t skipAmore=kFALSE;
146 
147  if ( config.GetConfigurationMap()->GetValue("LaserTriggerName") ) {
148  laserTriggerName=config.GetConfigurationMap()->GetValue("LaserTriggerName")->GetName();
149  printf("TPCCEda: Laser trigger class name set to: %s.\n",laserTriggerName.Data());
150  }
151 
152  if ( config.GetConfigurationMap()->GetValue("MonitoringType") ) {
153  monitoringType=config.GetConfigurationMap()->GetValue("MonitoringType")->GetName();
154  printf("TPCCEda: Monitoring type set to: %s.\n",monitoringType.Data());
155  }
156 
157  if ( config.GetConfigurationMap()->GetValue("ForceLaserTriggerId") ) {
158  forceTriggerId=TMath::Nint(config.GetValue("ForceLaserTriggerId"));
159  printf("TPCCEda: Only processing triggers with Id: %d.\n",forceTriggerId);
160  }
161 
162  if ( config.GetConfigurationMap()->GetValue("SaveOption") ) {
163  saveOption=TMath::Nint(config.GetValue("SaveOption"));
164  printf("TPCCEda: Saving option set to: %d.\n",saveOption);
165  }
166 
167  if ( config.GetConfigurationMap()->GetValue("SkipAmore") ) {
168  skipAmore=((TObjString*)config.GetConfigurationMap()->GetValue("SkipAmore"))->GetString().Atoi();
169  printf("TPCCEda: Skip Amore set in config\n");
170  }
171 
172  //subsribe to laser triggers only in physics partition
173  //if the trigger class is not available the return value is -1
174  //in this case we are most probably running as a standalone
175  // laser run and should request all events
176  unsigned char classIdptr=0;
177  int retClassId=daqDA_getClassIdFromName(laserTriggerName.Data(),&classIdptr);
178  if (retClassId==0){
179  //interleaved laser in physics runs
180  //select proper trigger class id
181  char c[5];
182  snprintf(c,sizeof(c),"%u",(unsigned int)classIdptr);
183  char *table[5] = {"PHY",const_cast<char*>(monitoringType.Data()),"*",c,NULL};
184  monitorDeclareTableExtended(table);
185  printf("TPCCEda: Using monitoring table: (PHY, %s, *, %s)\n",monitoringType.Data(),c);
186  } else if (retClassId==-1){
187  //global partition without laser triggered events
188  //the DA should exit properly without processing
189  printf("TPCCEda: Laser trigger class '%s' was not found among trigger class names. Will stop processing.\n",laserTriggerName.Data());
190  return 0;
191  } else if (retClassId==-2){
192  //standalone case, accept all physics events
193  char *table[5] = {"PHY","Y","*","*",NULL};
194  monitorDeclareTableExtended(table);
195  printf("TPCCEda: Using all trigger class Ids\n");
196  } else {
197  printf("TPCCEda: Unknown return value of 'daqDA_getClassIdFromName': %d\n",retClassId);
198  return -2;
199  }
200 
201  //see if we should force the trigger id
202  if (forceTriggerId>-1){
203  char c[5];
204  sprintf(c,"%d",forceTriggerId);
205  char *table[5] = {"PHY","Y","*",c,NULL};
206  monitorDeclareTableExtended(table);
207  }
208 
209 
210  // if test setup get parameters from $DAQDA_TEST_DIR
211  if (!mapping){
212  /* copy locally the mapping file from daq detector config db */
213  sprintf(localfile,"./%s",MAPPING_FILE);
214  status = daqDA_DB_getFile(MAPPING_FILE,localfile);
215  if (status) {
216  printf("TPCCEda: Failed to get mapping file (%s) from DAQdetDB, status=%d\n", MAPPING_FILE, status);
217  return -1;
218  }
219 
220  /* open the mapping file and retrieve mapping object */
221  TFile *fileMapping = new TFile(MAPPING_FILE, "read");
222  mapping = (AliTPCmapper*) fileMapping->Get("tpcMapping");
223  delete fileMapping;
224  }
225 
226  if (mapping == 0) {
227  printf("TPCCEda: Failed to get mapping object from %s. ...\n", MAPPING_FILE);
228  return -1;
229  } else {
230  printf("TPCCEda: Got mapping object from %s\n", MAPPING_FILE);
231  }
232 
233 
234  //create calibration object
235  AliTPCCalibCE *calibCE=new AliTPCCalibCE(config.GetConfigurationMap()); // central electrode calibration
236  calibCE->SetAltroMapping(mapping->GetAltroMapping()); // Use altro mapping we got from daqDetDb
237 
238  //amore update interval
239  Double_t updateInterval=300; //seconds
240  Double_t valConf=config.GetValue("AmoreUpdateInterval");
241  if ( valConf>0 ) updateInterval=valConf;
242  //timer
243  TStopwatch stopWatch;
244 
245  //===========================//
246  // loop over RAW data files //
247  //==========================//
248  int nevents=0;
249  int neventsOld=0;
250  size_t counter=0;
251  for ( i=1; i<argc; i++) {
252 
253  /* define data source : this is argument i */
254  printf("TPCCEda: Processing file %s\n", argv[i]);
255  status=monitorSetDataSource( argv[i] );
256  if (status!=0) {
257  printf("TPCCEda: monitorSetDataSource() failed : %s\n",monitorDecodeError(status));
258  return -1;
259  }
260 
261  Bool_t hasNewData=kFALSE;
262  /* read until EOF */
263  while (true) {
264  struct eventHeaderStruct *event;
265 
266  /* check shutdown condition */
267  if (daqDA_checkShutdown()) {break;}
268 
269  /* get next event (blocking call until timeout) */
270  status=monitorGetEventDynamic((void **)&event);
271  if (status==MON_ERR_EOF) {
272  printf ("TPCCEda: End of File %d detected\n",i);
273  break; /* end of monitoring file has been reached */
274  }
275 
276  if (status!=0) {
277  printf("TPCCEda: monitorGetEventDynamic() failed : %s\n",monitorDecodeError(status));
278  break;
279  }
280 
281  /* retry if got no event */
282  if (event==NULL){
283  //use time in between bursts to
284  // send the data to AMOREdb
285  if (stopWatch.RealTime()>updateInterval && hasNewData){
286  calibCE->Analyse();
287  if (!skipAmore) SendToAmoreDB(calibCE,runNb);
288  stopWatch.Start();
289  } else {
290  stopWatch.Continue();
291  }
292  //debug output
293  if (nevents>neventsOld){
294  printf ("TPCCEda: %d events processed, %d used\n",nevents,calibCE->GetNeventsProcessed());
295  neventsOld=nevents;
296  }
297  hasNewData=kFALSE;
298  continue;
299  }
300 
301  /* skip start/end of run events */
302  if ( (event->eventType != physicsEvent) && (event->eventType != calibrationEvent) ){
303  free(event);
304  continue;
305  }
306 
307 
308  // get the run number
309  runNb = event->eventRunNb;
310 
311  // CE calibration
312  calibCE->ProcessEvent(event);
313  hasNewData=kTRUE;
314 
315  /* free resources */
316  free(event);
317  ++nevents;
318  }
319  }
320 
321  //
322  // Analyse CE data and write them to rootfile
323  //
324  printf ("TPCCEda: %d events processed, %d used\n",nevents,calibCE->GetNeventsProcessed());
325 
326  //save data to file
327  calibCE->DumpToFile(RESULT_FILE,Form("name=tpcCalibCE,type=%d",saveOption));
328  printf("TPCCEda: Wrote %s\n",RESULT_FILE);
329 
330  /* store the result file on FES */
331  status=daqDA_FES_storeFile(RESULT_FILE,FILE_ID);
332  if (status) {
333  status = -2;
334  }
335 
336  if (!skipAmore){
337  printf("TPCCEda: Amore part\n");
338  calibCE->Analyse();
339  SendToAmoreDB(calibCE,runNb);
340  }
341 
342  delete calibCE;
343  return status;
344 }
345 
346 
347 void SendToAmoreDB(AliTPCCalibCE *calibCE, unsigned long32 runNb)
348 {
349  //AMORE
350 // printf ("AMORE part\n");
351  const char *amoreDANameorig=gSystem->Getenv("AMORE_DA_NAME");
352  //cheet a little -- temporary solution (hopefully)
353  //
354  //currently amoreDA uses the environment variable AMORE_DA_NAME to create the mysql
355  //table in which the calib objects are stored. This table is dropped each time AmoreDA
356  //is initialised. This of course makes a problem if we would like to store different
357  //calibration entries in the AMORE DB. Therefore in each DA which writes to the AMORE DB
358  //the AMORE_DA_NAME env variable is overwritten.
359  gSystem->Setenv("AMORE_DA_NAME",Form("TPC-%s",FILE_ID));
360  //
361  // end cheet
362  TGraph *grA=calibCE->MakeGraphTimeCE(-1,0,2);
363  TGraph *grC=calibCE->MakeGraphTimeCE(-2,0,2);
364  TDatime time;
365  TObjString info(Form("Run: %u; Date: %s",runNb,time.AsSQLString()));
366  amore::da::AmoreDA amoreDA(amore::da::AmoreDA::kSender);
367  Int_t statusDA=0;
368  statusDA+=amoreDA.Send("CET0",calibCE->GetCalPadT0());
369  statusDA+=amoreDA.Send("CEQ",calibCE->GetCalPadQ());
370  statusDA+=amoreDA.Send("CERMS",calibCE->GetCalPadRMS());
371  statusDA+=amoreDA.Send("DriftA",grA);
372  statusDA+=amoreDA.Send("DriftC",grC);
373  statusDA+=amoreDA.Send("Info",&info);
374  if ( statusDA!=0 )
375  printf("TPCCEda: Waring: Failed to write one of the calib objects to the AMORE database\n");
376  // reset env var
377  if (amoreDANameorig) gSystem->Setenv("AMORE_DA_NAME",amoreDANameorig);
378  if (grA) delete grA;
379  if (grC) delete grC;
380 }
381 
382 
printf("Chi2/npoints = %f\n", TMath::Sqrt(chi2/npoints))
const TObjArray * GetCalPadRMS() const
Definition: AliTPCCalibCE.h:57
#define CONFIG_FILE
Definition: TPCCEda.cxx:41
static void SetClassDebugLevel(const char *className, Int_t level)
Definition: AliLog.cxx:513
const TMap * GetConfigurationMap() const
AliTPCmapper.
Definition: AliTPCmapper.h:15
const TObjArray * GetCalPadT0() const
Definition: AliTPCCalibCE.h:54
TROOT * gROOT
#define MAPPING_FILE
Definition: TPCCEda.cxx:40
Bool_t ProcessEvent(AliTPCRawStreamV3 *const rawStreamV3)
virtual void DumpToFile(const Char_t *filename, const Char_t *dir="", Bool_t append=kFALSE)
virtual void Analyse()
AliExternalInfo info
Float_t GetValue(const char *name) const
Int_t GetNeventsProcessed() const
static void SetModuleDebugLevel(const char *module, Int_t level)
Definition: AliLog.cxx:485
Class for Parsing simple text configuration files.
const TObjArray * GetCalPadQ() const
Definition: AliTPCCalibCE.h:56
Implementation of the TPC Central Electrode calibration.
Definition: AliTPCCalibCE.h:29
#define FILE_ID
Definition: TPCCEda.cxx:39
AliTPCAltroMapping ** GetAltroMapping()
Definition: AliTPCmapper.h:28
void SetAltroMapping(AliTPCAltroMapping **mapp)
int main(int argc, char **argv)
Definition: TPCCEda.cxx:91
TGraph * MakeGraphTimeCE(Int_t sector, Int_t xVariable=0, Int_t fitType=0, Int_t fitParameter=0)
void SendToAmoreDB(AliTPCCalibCE *calibCE, unsigned long32 runNb)
Definition: TPCCEda.cxx:347
#define RESULT_FILE
Definition: TPCCEda.cxx:38