2 /**************************************************************************
3 * This file is property of and copyright by the ALICE HLT Project *
4 * ALICE Experiment at CERN, All rights reserved. *
6 * Primary Authors: Jochen Thaeder <thaeder@kip.uni-heidelberg.de> *
7 * for The ALICE HLT Project. *
9 * Permission to use, copy, modify and distribute this software and its *
10 * documentation strictly for non-commercial purposes is hereby granted *
11 * without fee, provided that the above copyright notice appears in all *
12 * copies and that both the copyright notice and this permission notice *
13 * appear in the supporting documentation. The authors make no claims *
14 * about the suitability of this software for any purpose. It is *
15 * provided "as is" without express or implied warranty. *
16 **************************************************************************/
18 /** @file AliHLTTPCCalibTimeComponent.cxx
19 @author Kalliopi Kanaki
21 @brief A calibration component for interfacing the offline calculation of TPC drift velocity correction
24 // see header file for class documentation
26 // refer to README to build package
28 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
35 #include "AliHLTTPCCalibTimeComponent.h"
36 #include "AliHLTTPCDefinitions.h"
37 #include "AliHLTMisc.h"
39 #include "AliESDEvent.h"
40 #include "AliESDtrack.h"
41 #include "AliESDfriend.h"
43 #include "AliTPCcalibTime.h"
44 #include "AliTPCcalibCalib.h"
45 #include "AliTPCseed.h"
46 #include "AliTPCcalibDB.h"
47 #include "AliTPCClusterParam.h"
49 #include "AliHLTTPCOfflineCluster.h"
50 #include "AliHLTTPCSpacePointData.h"
51 #include "AliHLTTPCTrackletDataFormat.h"
52 #include "AliHLTExternalTrackParam.h"
53 #include "AliHLTGlobalBarrelTrack.h"
54 #include "AliHLTTPCTransform.h"
56 #include "TObjArray.h"
60 #include "THnSparse.h"
61 #include "TGraphErrors.h"
66 #include "AliHLTReadoutList.h"
68 ClassImp(AliHLTTPCCalibTimeComponent) // ROOT macro for the implementation of ROOT specific class methods
70 AliHLTTPCCalibTimeComponent::AliHLTTPCCalibTimeComponent()
80 // see header file for class documentation
82 // refer to README to build package
84 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
86 for(int i=0; i<fkNPartition; i++){
87 fPartitionClusters[i] = 0;
88 fNPartitionClusters[i] = 0;
92 const char* AliHLTTPCCalibTimeComponent::fgkOCDBEntry="HLT/ConfigTPC/TPCCalibTime";
94 AliHLTTPCCalibTimeComponent::~AliHLTTPCCalibTimeComponent(){
95 // see header file for class documentation
97 for(int i=0; i<fkNPartition; i++){
98 delete[] fPartitionClusters[i];
102 const char* AliHLTTPCCalibTimeComponent::GetComponentID() {
103 // see header file for class documentation
104 return "TPCCalibTime";
107 void AliHLTTPCCalibTimeComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list){
108 // see header file for class documentation
111 list.push_back( kAliHLTDataTypeTObjArray|kAliHLTDataOriginTPC ); // output of TPCCalibSeedMaker
112 list.push_back( AliHLTTPCDefinitions::fgkClustersDataType ); // output of the TPC CF
113 list.push_back( kAliHLTDataTypeTrack|kAliHLTDataOriginTPC ); // output of the global merger
114 list.push_back( kAliHLTDataTypeESDObject|kAliHLTDataOriginOut ); // output of global esd converter
117 AliHLTComponentDataType AliHLTTPCCalibTimeComponent::GetOutputDataType() {
118 // see header file for class documentation
120 return AliHLTTPCDefinitions::CalibCEDataType()|kAliHLTDataOriginOut;
123 void AliHLTTPCCalibTimeComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier ){
124 // see header file for class documentation
126 constBase = fOutputSize;
127 inputMultiplier = 0; // to be estimated
130 AliHLTComponent* AliHLTTPCCalibTimeComponent::Spawn(){
131 // see header file for class documentation
133 return new AliHLTTPCCalibTimeComponent();
137 Int_t AliHLTTPCCalibTimeComponent::ScanConfigurationArgument( Int_t argc, const char** argv ){
138 // see header file for class documentation
140 if (argc<=0) return 0;
142 TString argument=argv[i];
145 if (argument.CompareTo("-output-size")==0) {
146 if (++i>=argc) return -EPROTO;
148 fOutputSize=argument.Atoi();
154 Int_t AliHLTTPCCalibTimeComponent::InitCalibration(){
155 // see header file for class documentation
157 AliTPCcalibDB::Instance()->SetRun(AliHLTMisc::Instance().GetCDBRunNo());
158 AliTPCcalibDB::Instance()->GetClusterParam()->SetInstance(AliTPCcalibDB::Instance()->GetClusterParam());
161 // AliTPCcalibDB *calib = AliTPCcalibDB::Instance();
164 // HLTError("AliTPCcalibDB does not exist");
168 // AliTPCClusterParam *clusPar = calib->GetClusterParam();
170 // HLTError("OCDB entry TPC/Calib/ClusterParam (AliTPCcalibDB::GetClusterParam()) is not available.");
174 // first configure the default
176 if (iResult>=0) iResult=ConfigureFromCDBTObjString(fgkOCDBEntry);
178 if(fCalibTime) return EINPROGRESS;
179 fCal = new AliTPCcalibCalib();
181 fSeedArray = new TObjArray();
187 Int_t AliHLTTPCCalibTimeComponent::DeinitCalibration() {
188 // see header file for class documentation
190 if(fCalibTime) delete fCalibTime; fCalibTime = NULL;
191 if(fCal) delete fCal; fCal = NULL;
192 if(fSeedArray) delete fSeedArray; fSeedArray = NULL;
194 //if(fESDfriend) delete fESDfriend; fESDfriend = NULL;
196 //if(arr) delete arr; arr = NULL;
201 int AliHLTTPCCalibTimeComponent::Reconfigure(const char* cdbEntry, const char* /*chainId*/){
202 // see header file for class documentation
204 // configure from the specified antry or the default one
205 const char* entry=cdbEntry;
206 if(!entry || entry[0]==0){
210 return ConfigureFromCDBTObjString(entry);
213 Int_t AliHLTTPCCalibTimeComponent::ProcessCalibration( const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/ ){
214 // see header file for class documentation
216 if(GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR )) return 0;
218 const AliHLTComponentBlockData *iter = NULL;
220 //--------------- output over TObjArray of AliTPCseed objects (output of TPCSeedMaker) -------------------//
222 // A previous component in the chain (TPCSeedMaker) has processed the TPC clusters and tracks and created a TClonesArray of AliTPCseed objects
223 // In this loop the iterator accesses this array stored in memory, in order to use it in the next loop over the AliESDevent of the HLT
225 // for(iter = (TObject*)GetFirstInputObject(kAliHLTDataTypeTObjArray|kAliHLTDataOriginTPC); iter != NULL; iter = (TObject*)GetNextInputObject()){
227 // if(GetDataType(iter) != (kAliHLTDataTypeTObjArray | kAliHLTDataOriginTPC)) continue;
228 // fSeedArray = dynamic_cast<TClonesArray*>(iter);
231 // int nInputClusters = 0;
232 // int nInputTracks = 0;
234 //TObjArray *arr = new TObjArray(1000);
235 //arr->SetOwner(kTRUE);
239 for( int i=0; i<fkNPartition; i++ ){
240 delete[] fPartitionClusters[i];
241 fPartitionClusters[i] = 0;
242 fNPartitionClusters[i] = 0;
247 //------------------- loop over clusters -------------//
249 for(iter = GetFirstInputBlock(AliHLTTPCDefinitions::fgkClustersDataType); iter != NULL; iter = GetNextInputBlock()){
251 if ( iter->fDataType != AliHLTTPCDefinitions::fgkClustersDataType ) continue;
253 Int_t slice = AliHLTTPCDefinitions::GetMinSliceNr(iter->fSpecification);
254 Int_t partition = AliHLTTPCDefinitions::GetMinPatchNr(iter->fSpecification);
256 Int_t slicepartition = slice*6+partition;
258 if(slicepartition > fkNPartition){
259 HLTWarning("Wrong header of TPC cluster data, slice %d, partition %d", slice, partition );
263 AliHLTTPCClusterData* inPtrSP = ( AliHLTTPCClusterData* )( iter->fPtr );
264 // nInputClusters += inPtrSP->fSpacePointCnt;
266 delete[] fPartitionClusters[slicepartition];
267 fPartitionClusters[slicepartition] = new AliTPCclusterMI[inPtrSP->fSpacePointCnt];
268 fNPartitionClusters[slicepartition] = inPtrSP->fSpacePointCnt;
270 // create offline clusters out of the HLT clusters
271 // todo: check which cluster information is really needed for the dEdx
272 for(unsigned int i = 0; i < inPtrSP->fSpacePointCnt; i++){
273 AliHLTTPCSpacePointData *chlt = &( inPtrSP->fSpacePoints[i] );
274 AliTPCclusterMI *c = fPartitionClusters[slicepartition]+i;
278 c->SetSigmaY2(chlt->fSigmaY2);
280 c->SetSigmaZ2(chlt->fSigmaZ2);
281 c->SetQ( chlt->fCharge );
282 c->SetMax( chlt->fQMax );
284 Float_t padtime[3] = {0,chlt->fY,chlt->fZ};
285 AliHLTTPCTransform::Slice2Sector(slice,chlt->fPadRow, sector, row);
286 AliHLTTPCTransform::Local2Raw( padtime, sector, row);
287 c->SetDetector( sector );
289 c->SetPad( (Int_t) padtime[1] );
290 c->SetTimeBin( (Int_t) padtime[2] );
292 } // end of loop over blocks of clusters
297 //---------- loop over merged tracks ------------------ //
299 for(const AliHLTComponentBlockData *pBlock = GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC); pBlock != NULL; pBlock = GetNextInputBlock()){
301 AliHLTTracksData *dataPtr = (AliHLTTracksData*) pBlock->fPtr;
302 nTracks = dataPtr->fCount;
304 AliHLTExternalTrackParam *currTrack = dataPtr->fTracklets;
306 //nInputTracks += nTracks;
308 for(int itr=0; itr<nTracks && ( (AliHLTUInt8_t *)currTrack < ((AliHLTUInt8_t *) pBlock->fPtr)+pBlock->fSize); itr++){
310 // create an offline track
311 AliHLTGlobalBarrelTrack gb(*currTrack);
313 tTPC.Set(gb.GetX(), gb.GetAlpha(), gb.GetParameter(), gb.GetCovariance());
315 // set the cluster pointers
316 for(UInt_t ic=0; ic<currTrack->fNPoints; ic++){
318 tTPC.SetNumberOfClusters(currTrack->fNPoints);
320 UInt_t id = currTrack->fPointIDs[ic];
321 int iSlice = AliHLTTPCSpacePointData::GetSlice(id);
322 int iPartition = AliHLTTPCSpacePointData::GetPatch(id);
323 int iCluster = AliHLTTPCSpacePointData::GetNumber(id);
325 if(iSlice<0 || iSlice>36 || iPartition<0 || iPartition>5){
326 HLTError("Corrupted TPC cluster Id: slice %d, partition %d, cluster %d", iSlice, iPartition, iCluster);
330 AliTPCclusterMI *partitionClusters = fPartitionClusters[iSlice*6 + iPartition];
332 if(!partitionClusters){
333 HLTError("Clusters are missed for slice %d, partition %d", iSlice, iPartition );
337 if(iCluster >= fNPartitionClusters[iSlice*6 + iPartition]){
338 HLTError("TPC slice %d, partition %d: ClusterID==%d >= N Clusters==%d ", iSlice, iPartition,iCluster, fNPartitionClusters[iSlice*6 + iPartition]);
342 AliTPCclusterMI *c = &(partitionClusters[iCluster]);
343 int sec = c->GetDetector();
344 int row = c->GetRow();
345 if(sec >= 36) row = row + AliHLTTPCTransform::GetNRowLow();
347 tTPC.SetClusterPointer(row, c);
349 AliTPCTrackerPoint &point = *( tTPC.GetTrackPoint( row ) );
350 //tTPC.Propagate( TMath::DegToRad()*(sec%18*20.+10.), c->GetX(), fSolenoidBz );
351 Double_t angle2 = tTPC.GetSnp()*tTPC.GetSnp();
352 angle2 = (angle2<1) ?TMath::Sqrt(angle2/(1-angle2)) :10.;
353 point.SetAngleY( angle2 );
354 point.SetAngleZ( tTPC.GetTgl() );
355 } // end of associated cluster loop
357 AliTPCseed *seed = &(tTPC);
358 fSeedArray->Add(seed);
360 unsigned int step = sizeof( AliHLTExternalTrackParam ) + currTrack->fNPoints * sizeof( unsigned int );
361 currTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currTrack) + step );
363 }// end of vector track loop
364 } // end of loop over blocks of merged tracks
366 HLTInfo("Number of reconstructed tracks %d, number of produced seeds %d\n", nTracks, fSeedArray->GetEntries());
369 //----------- loop over output of global esd converter ----------------//
371 // In this loop we access the AliESDevent that was produced by the HLT and is stored in memory. There should exist 1 object
372 // of type kAliHLTDataTypeESDObject per event.
374 TObject *iterOb = NULL;
375 for(iterOb = (TObject*)GetFirstInputObject(kAliHLTDataTypeESDObject | kAliHLTDataOriginOut); iterOb != NULL; iterOb = (TObject*)GetNextInputObject()){
377 if(GetDataType(iterOb) != (kAliHLTDataTypeESDObject | kAliHLTDataOriginOut)) continue;
379 fESDevent = dynamic_cast<AliESDEvent*>(iterOb);
380 fESDevent->GetStdContent();
382 HLTInfo("Number of seeds: %i\n", fSeedArray->GetEntriesFast()); // access of the info from the previous loop over the AliTPCseed array
384 fCal->UpdateEventInfo(fESDevent);
385 for(Int_t i=0; i<fSeedArray->GetEntriesFast(); i++){ // loop over TObjArray with seeds
387 AliTPCseed *seed = (AliTPCseed*)fSeedArray->UncheckedAt(i);
388 fESDtrack = fESDevent->GetTrack(i);
389 if(!fESDtrack || !seed) continue;
391 //if(fESDtrack->GetID() != seed->GetLabel()) {
392 // HLTWarning("Mismatch of track id between seed and ESD track: %i, %i\n", fESDtrack->GetID(), seed->GetLabel());
396 if(seed->GetNumberOfClusters()==0) continue;
397 fCal->RefitTrack(fESDtrack, seed, GetBz()); // update AliESDtrack and AliTPCseed info, acccording to Marian's request
399 AliTPCseed *seedCopy = new AliTPCseed(*seed, kTRUE);
400 fESDtrack->AddCalibObject(seedCopy); // add the AliTPCseed as a friend track to the AliESDtrack (to be accessed in TPC/AliTPCcalibTime.cxx)
402 //fESDfriendTrack = const_cast<AliESDfriendTrack*>(fESDtrack->GetFriendTrack());
406 Int_t startTime = fESDevent->GetTimeStamp()-60*60*1; //Start time one hour before first event, will make precise cuts later.
407 Int_t endTime = fESDevent->GetTimeStamp()+60*60*23; //End time 23 hours after first event.
408 // create the calibration object that will call the offline functions
409 fCalibTime = new AliTPCcalibTime("calibTime","time dependent Vdrift calibration", startTime, endTime, 20*60);
411 printf("fCalibTime pointer is NULL\n");
415 fCalibTime->SetStreamLevel(20);
416 fCalibTime->SetDebugLevel(20);
417 printf("fCalibTime = %i, startTime = %i, endTime = %i \n", fCalibTime!=0, startTime, endTime);
418 fESDfriend = new AliESDfriend();
419 fESDevent->GetESDfriend(fESDfriend);
420 fESDevent->SetESDfriend(fESDfriend);
421 fESDevent->AddObject(fESDfriend);
422 // create the AliESDfriend and add it to the event, now both the friend tracks and the friends are available for the offline functions to be called
424 fCalibTime->UpdateEventInfo(fESDevent); // needed for getting the run number and time stamp information correct on the offline side
425 fCalibTime->Process(fESDevent); // first offline function called
427 // delete fESDfriend;
429 //PushBack( (TObject*)fCalibTime, AliHLTTPCDefinitions::fgkCalibCEDataType | kAliHLTDataOriginOut, 0x0);
434 Int_t AliHLTTPCCalibTimeComponent::ShipDataToFXS( const AliHLTComponentEventData& /*evtData*/, AliHLTComponentTriggerData& /*trigData*/ ){
435 // see header file for class documentation
437 HLTInfo("Shipping data to FXS...\n");
439 fCalibTime->Analyze(); // called at the end of the run or event modulo
441 // the rest of the histogram and graph declarations were copied by Dag as a first attempt to get the start/end time bin "automatically". Perhaps we need some more
442 // thinking here to avoid copying all these lines that might chain in offline without HLT realizing.
444 THnSparse* addHist = fCalibTime->GetHistoDrift("all");
445 if(!addHist) return -1;
447 //Identifying used range of histogram
449 Int_t startTimeBin = 0;
450 Int_t endTimeBin = 0;
452 TH1D *histoTime = addHist->Projection(0);
454 startTimeBin = histoTime->FindFirstBinAbove(0);
455 endTimeBin = histoTime->FindLastBinAbove(0);
456 printf("startTimeBin = %i endTimeBin = %i\n", startTimeBin, endTimeBin);
457 printf("startTimeBinCentre = %f endTimeBinCentre = %f\n", histoTime->GetBinCenter(startTimeBin), histoTime->GetBinCenter(endTimeBin));
458 printf("startTimeBinWidth = %f endTimeBinWidth = %f\n", histoTime->GetBinWidth(startTimeBin), histoTime->GetBinWidth(endTimeBin));
459 delete histoTime; histoTime = 0;
462 Int_t startPtBin = 0;
464 TH1D *histoPt = addHist->Projection(1);
466 startPtBin = histoPt->FindFirstBinAbove(0);
467 endPtBin = histoPt->FindLastBinAbove(0);
468 printf("startPtBin = %i endPtBin = %i\n", startPtBin, endPtBin);
469 printf("startPtBinCentre = %f endPtBinCentre = %f\n", histoPt->GetBinCenter(startPtBin), histoPt->GetBinCenter(endPtBin));
470 printf("startPtinWidth = %f endPtBinWidth = %f\n", histoPt->GetBinWidth(startPtBin), histoPt->GetBinWidth(endPtBin));
471 delete histoPt; histoPt = 0;
474 Int_t startVdBin = 0;
476 TH1D *histoVd = addHist->Projection(2);
478 startVdBin = histoVd->FindFirstBinAbove(0);
479 endVdBin = histoVd->FindLastBinAbove(0);
480 printf("startVdBin = %i endVdBin = %i\n", startVdBin, endVdBin);
481 printf("startVdBinCentre = %f endVdBinCentre = %f\n", histoVd->GetBinCenter(startVdBin), histoVd->GetBinCenter(endVdBin));
482 printf("startVdBinWidth = %f endVdBinWidth = %f\n", histoVd->GetBinWidth(startVdBin), histoVd->GetBinWidth(endVdBin));
483 delete histoVd; histoVd = 0;
486 Int_t startRunBin = 0;
488 TH1D *histoRun = addHist->Projection(3);
490 startRunBin = histoRun->FindFirstBinAbove(0);
491 endRunBin = histoRun->FindLastBinAbove(0);
492 printf("startRunBin = %i endRunBin = %i\n", startRunBin, endRunBin);
493 printf("startRunBinCentre = %f endRunBinCentre = %f\n", histoRun->GetBinCenter(startRunBin), histoRun->GetBinCenter(endRunBin));
494 printf("startRunBinWidth = %f endRunBinWidth = %f\n", histoRun->GetBinWidth(startRunBin), histoRun->GetBinWidth(endRunBin));
495 delete histoRun; histoRun = 0;
498 TObjArray *vdriftArray = new TObjArray();
499 if(!vdriftArray) return -2;
501 TObjArray *array = fCalibTime->GetHistoDrift();
502 if(!array) return -3;
504 TIterator *iterator = array->MakeIterator();
505 if(!iterator) return -4;
508 THnSparse *hist = NULL;
509 while((hist = (THnSparseF*)iterator->Next())){
511 //if(!hist) continue;
513 hist->GetAxis(0)->SetRange(startTimeBin, endTimeBin);
514 hist->GetAxis(1)->SetRange(startPtBin, endPtBin);
515 hist->GetAxis(0)->SetRange(startVdBin, endVdBin);
516 hist->GetAxis(3)->SetRange(startRunBin, endRunBin);
518 TString name = hist->GetName();
519 Int_t dim[4] = {0,1,2,3};
520 THnSparse *newHist = hist->Projection(4,dim);
521 newHist->SetName(name);
522 vdriftArray->Add(newHist);
524 TGraphErrors *graph = AliTPCcalibBase::FitSlices(newHist,2,0,400,100,0.05,0.95, kTRUE);
525 printf("name = %s graph = %i\n", name.Data(), graph==0);
526 if(!graph || !graph->GetN()) continue;
527 printf("name = %s graph = %i, N = %i\n", name.Data(), graph==0, graph->GetN());
528 Int_t pos = name.Index("_");
529 name = name(pos,name.Capacity()-pos);
530 TString graphName = graph->ClassName();
533 graph->SetName(graphName);
534 printf("name = %s\n", graphName.Data());
535 vdriftArray->Add(graph);
537 //Currently, AliSplineFits can not be given names...
538 //AliSplineFit* fit=new AliSplineFit();
539 //fit->SetGraph(graph);
540 //fit->SetMinPoints(graph->GetN()+1);
541 //fit->InitKnots(graph,2,0,0.001);
543 //TString fiName=fit->ClassName();
547 //fit->SetName(fiName.Data());
548 //printf("name=%s\n", fiName.Data());
549 //vdriftArray->Add(fit);
552 THnSparse *laserHist = NULL;
553 TGraphErrors *laserGraph = NULL;
554 TString laserName = "";
556 //Histograms and graphs for A side lasers
557 laserHist = fCalibTime->GetHistVdriftLaserA(1);
560 laserName=laserHist->ClassName();
561 laserName+="_MEAN_DRIFT_LASER_ALL_A";
563 laserHist->SetName(laserName);
564 vdriftArray->Add(laserHist);
565 laserGraph=AliTPCcalibBase::FitSlices(laserHist,2,0,400,100,0.05,0.95, kTRUE);
566 if(laserGraph && laserGraph->GetN()){
567 laserName=laserGraph->GetName();
568 laserName+="_MEAN_DRIFT_LASER_ALL_A";
570 laserGraph->SetName(laserName);
571 vdriftArray->Add(laserGraph);
575 //Histograms and graphs for C side lasers
576 laserHist=fCalibTime->GetHistVdriftLaserC(1);
578 laserName=laserHist->ClassName();
579 laserName+="_MEAN_DRIFT_LASER_ALL_C";
581 laserHist->SetName(laserName);
582 vdriftArray->Add(laserHist);
583 laserGraph=AliTPCcalibBase::FitSlices(laserHist,2,0,400,100,0.05,0.95, kTRUE);
584 if(laserGraph && laserGraph->GetN()){
585 laserName=laserGraph->GetName();
586 laserName+="_MEAN_DRIFT_LASER_ALL_C";
588 laserGraph->SetName(laserName);
589 vdriftArray->Add(laserGraph);
593 //Meatdata set in off-line...
594 //AliCDBMetaData *metaData= new AliCDBMetaData();
595 //metaData->SetObjectClassName("TObjArray");
596 //metaData->SetResponsible("Dag Toppe Larsen");
597 //metaData->SetBeamPeriod(1);
598 //metaData->SetAliRootVersion("05-25-01"); //root version
599 //metaData->SetComment("Calibration of the time dependence of the drift velocity due to pressure and temperature changes");
600 //AliCDBId* id1=NULL;
601 //if(end) id1=new AliCDBId("TPC/Calib/TimeDrift", runNumber, end);
602 //else id1=new AliCDBId("TPC/Calib/TimeDrift", runNumber, runNumber);
603 //AliCDBStorage* gStorage = AliCDBManager::Instance()->GetStorage("local://$ALICE_ROOT/OCDB");
604 //gStorage->Put(vdriftArray, (*id1), metaData);
605 //printf("done runNumber=%i, end=%i\n", runNumber, end);
607 static AliHLTReadoutList rdList(AliHLTReadoutList::kTPC);
610 TFile *file = TFile::Open("vdrift.root", "RECREATE");
611 vdriftArray->Write();
615 file = TFile::Open("calibTime.root", "RECREATE");
620 // the vdriftArray is pushed to the HLT-FXSsubscriber
621 PushToFXS( (TObject*)vdriftArray, "TPC", "TIMEDRIFT", &rdList );
623 //Should array be deleted now?
625 // vdriftArray.Clear();
626 // delete vdriftArray;