8abfdb64ea49ea924f46385dcfef924c6b18cb3e
[u/mrichter/AliRoot.git] / STEER / CDB / AliOCDBtoolkit.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /*
17
18   Primary goal of the proposal was to provide functionality to browse and compare the content of the OCDB
19   specified by different means.
20
21   a.) galice.root               - as currently implemented by Ruben in MC (or any file with cdbMap, and cdbList)
22   b.) AliESDs.root              - for the reconstructed data
23   c.) ocdb snapshot             - as used for grid productions
24   d.) TMap(s)                   - as used internally in galice.root and AliESDs,root  
25   e.) log file (if possible)    - looks logs aways used similar syntax, tested and working
26   f.) C macro                   - custom macro
27
28   Content comparison should be done:
29   a.) on the level of symbolic links 
30   b.) on the level of content itself 
31   - by by byte comparison dif
32   - data member by data member comparison
33
34   Implementation assumption:
35   All input formats (a .. f) will  be converted to the TMap storages and TList if AliCDBIds 
36
37   Example usage:
38   AliOCDBtoolkit::MakeDiffExampleUseCase();
39   or from the AliOCDBtoolkit.sh in propmpt
40   ocdbMakeTable AliESDs.root ESD OCDBrec.list
41   ocdbMakeTable galice.root MC OCDBsim.list
42
43   
44    
45
46
47
48   //=============================================================================
49   // Functionality to dump content of  objects in human readable format
50   //=============================================================================
51   Use case examples 
52   1.) compare oontent of alignent OCDB files for differnt yers
53   2.) compare ClusterParam for different periods
54   
55   
56    
57   =================================================================================================================
58   // 1.)
59   // Compare alignment example:
60   // Compare TPC alignemnt 2013 and 2010
61   //
62   AliOCDBtoolkit::DumpOCDBFile("/cvmfs/alice.gsi.de/alice/data/2013/OCDB/TPC/Align/Data/Run0_999999999_v1_s0.root","TPCalign2013.dump",1,1);
63   AliOCDBtoolkit::DumpOCDBFile("/cvmfs/alice.gsi.de/alice/data/2010/OCDB/TPC/Align/Data/Run0_999999999_v1_s0.root","TPCalign2010.dump",1,1);
64   diff  TPCalign2013.dump TPCalign2010.dump > TPCalign2013_TPCalign2010.diff
65   //
66   //    
67   =================================================================================================================
68   //  2.) 
69   // Compare CluterParam OCDB etry
70   //
71   AliOCDBtoolkit::DumpOCDBFile("/cvmfs/alice.gsi.de/alice/data/2010/OCDB/TPC/Calib/ClusterParam/Run131541_999999999_v2_s0.root","2010_TPC_Calib_ClusterParam_Run131541_999999999_v2_s0.dump",1);
72   AliOCDBtoolkit:: AliOCDBtoolkit::DumpOCDBFile("/cvmfs/alice.gsi.de/alice/data/2010/OCDB/TPC/Calib/ClusterParam/Run0_999999999_v1_s0.root","2010_TPC_Calib_ClusterParam_Run0_999999999_v1_s0.dump",1);
73   AliOCDBtoolkit::DumpOCDBFile("/cvmfs/alice.gsi.de/alice/data/2013/OCDB/TPC/Calib/ClusterParam/Run0_999999999_v1_s0.root","2013_TPC_Calib_ClusterParam_Run0_999999999_v1_s0.dump",1);
74   diff 2010_TPC_Calib_ClusterParam_Run131541_999999999_v2_s0.dump 2010_TPC_Calib_ClusterParam_Run0_999999999_v1_s0.dump
75  
76
77 */
78
79 /*
80   To check:
81   1.) Verify hash value uasge as and MD5 sum - 
82
83  */
84
85 using namespace std;
86
87 // STD
88 #include <iostream>
89 #include <algorithm>
90 #include <sstream>
91 #include <stdexcept>
92 #include <functional>
93 #include "TRealData.h"
94 #include "TDataMember.h"
95 #include "TClass.h"
96 #include "TROOT.h"
97 #include <TVectorD.h>
98 //
99 #include "TSystem.h"
100 #include "TObjArray.h"
101 #include "TString.h"
102 #include "TTree.h"
103 #include "TMessage.h"
104 //
105 #include "AliCDBManager.h"
106 #include "AliCDBEntry.h"
107 #include "AliOCDBtoolkit.h"
108 #include "AliCDBStorage.h"
109
110
111 void AliOCDBtoolkit::MakeDiffExampleUseCase(){
112   //
113   // Example usage for the MC 
114   // To run example case, assuming presence of following files in working directory: 
115   //    - rec.log        
116   //    - galice.root   
117   //    - AliESDs.root
118   //
119   AliCDBManager * man = AliCDBManager::Instance();
120   AliOCDBtoolkit::LoadOCDBFromLog("rec.log",0);
121   const TMap *cdbMapLog= man->GetStorageMap();        // this is map of 
122   const TList *cdbListLog=man->GetRetrievedIds();     // this is list of AliCDBId
123   //  TList *cdbListLog0=man->GetRetrievedIds();     // this is list of AliCDBId
124   //
125   TFile *fmc = TFile::Open("galice.root");
126   TMap *cdbMapMC= (TMap*)fmc->Get("cdbMap");          // 
127   TList *cdbListMC0= (TList*)fmc->Get("cdbList");     // this is list of TObjStrings
128   TList *cdbListMC = AliOCDBtoolkit::ConvertListStringToCDBId(cdbListMC0);        // convert to the TObjArray of AliCDBids
129   //
130   TFile *fesd = TFile::Open("AliESDs.root");
131   TList *listESD = ((TTree*)fesd->Get("esdTree"))->GetUserInfo();
132   TMap *cdbMapESD= (TMap*)listESD->FindObject("cdbMap");  
133   TList *cdbListESD0= (TList*)listESD->FindObject("cdbList"); // this is list of TObjStrings
134   TList *cdbListESD = AliOCDBtoolkit::ConvertListStringToCDBId(cdbListESD0);              // convert to the TObjArray  of AliCDBids
135   //
136   //
137   //
138   printf("\n\n");
139   printf("Diff log>>>ESD\n\n:");
140   MakeDiff(cdbMapLog, cdbListLog, cdbMapESD, cdbListESD,0);
141   printf("\n\n");
142   printf("Diff ESD>>>log\n\n:");
143   MakeDiff(cdbMapESD, cdbListESD,cdbMapLog, cdbListLog,0);
144   // 
145   printf("\n\n");
146   printf("Diff ESD>>>MC\n\n:");
147   MakeDiff(cdbMapMC, cdbListMC, cdbMapESD, cdbListESD,0);
148 }
149
150
151 void AliOCDBtoolkit::DumpOCDBAsTxt(const TString fInput, const TString fType, const TString outfile){
152   //
153   //
154   //
155   TFile *file;
156   const TMap *cdbMap=0;
157   const TList *cdbList=0;
158   //
159   //
160   AliCDBManager * man = AliCDBManager::Instance();
161
162   if(fType.EqualTo("MC",TString::kIgnoreCase)){
163         file = TFile::Open(fInput.Data());
164         cdbMap = (TMap*)file->Get("cdbMap");
165         if (!cdbMap){
166           printf("cdbMap does not exist in input file\t%s. Exiting\n",fInput.Data());
167           return;
168         }
169         // 
170         man->SetDefaultStorage(((TPair*)cdbMap->FindObject("default"))->Value()->GetName());
171         TList *cdbListMC0 = (TList*)file->Get("cdbList");     // this is list of TObjStrings
172         cdbList = AliOCDBtoolkit::ConvertListStringToCDBId(cdbListMC0);        // convert to the TObjArray of AliCDBids
173   } 
174     else if(fType.EqualTo("ESD",TString::kIgnoreCase)){
175       file = TFile::Open(fInput.Data());
176       if (!file) {
177         printf("Input file  does not exist %s. Exiting\n",fInput.Data());
178         return;
179       }
180       TList *listESD = ((TTree*)file->Get("esdTree"))->GetUserInfo();
181       cdbMap = (TMap*)listESD->FindObject("cdbMap");  
182       if (!cdbMap){
183         printf("cdbMap does not exist in input file\t%s. Exiting\n",fInput.Data());
184         return;
185       }
186       man->SetDefaultStorage(((TPair*)cdbMap->FindObject("default"))->Value()->GetName());
187       TList *cdbListESD0= (TList*)listESD->FindObject("cdbList"); // this is list of TObjStrings
188       cdbList = ConvertListStringToCDBId(cdbListESD0);              // convert to the TObjArray  of AliCDBids
189     }
190     else if(fType.EqualTo("log",TString::kIgnoreCase)){
191         LoadOCDBFromLog(fInput.Data(),0);
192         cdbMap = man->GetStorageMap();        // this is map of 
193         cdbList =man->GetRetrievedIds();     // this is list of AliCDBId
194     }
195     else{
196         printf("unsupported option: %s",fType.Data());
197         return;
198     }
199   cout <<"BEGINDUMP:" << endl;
200   DumpOCDB(cdbMap,cdbList,outfile);
201 }
202
203
204 Bool_t AliOCDBtoolkit::ParseInfoFromOcdbString(TString ocdbString, TString &ocdbPath, Int_t &run0, Int_t &run1, Int_t &version, Int_t &subVersion){
205   // Functionalit
206   // Parse OCDB id string and provide basic ocdb information
207   //
208   //  a.) parse ocdbPath
209   Int_t indexBeginPath= ocdbString.Index("path: ")+7;
210   if (indexBeginPath<0) return kFALSE;
211   Int_t indexEndPath=ocdbString.Index(";",indexBeginPath);
212   if (indexEndPath<0) return kFALSE;
213   ocdbPath=TString(&(ocdbString.Data()[indexBeginPath]), indexEndPath-indexBeginPath-1);
214   // b.) parse runRange
215   Int_t indexRun0= ocdbString.Index(": [",indexEndPath)+3;
216   if (indexRun0<0) return kFALSE;
217   Int_t indexRun1= ocdbString.Index(",",indexRun0)+1;
218   if (indexRun1<0) return kFALSE;
219   run0=atoi(&(ocdbString.Data()[indexRun0]));
220   run1=atoi(&(ocdbString.Data()[indexRun1]));
221   AliCDBRunRange runRange(run0,run1);
222   //c.) parse version, subversion
223   Int_t indexVersion= ocdbString.Index("version: v",indexRun1)+10;
224   if (indexVersion<0) return kFALSE;
225   Int_t indexSubVersion= ocdbString.Index("_s",indexVersion)+2;
226   if (indexSubVersion<0) return kFALSE;
227   version=atoi(&(ocdbString.Data()[indexVersion]));
228   subVersion=atoi(&(ocdbString.Data()[indexSubVersion]));
229   return kTRUE;
230 }
231
232 Bool_t AliOCDBtoolkit::ParseInfoFromOcdbString(TString ocdbString, AliCDBId &cdbId){
233   //
234   // Parse OCDB id string and provide basic ocdb information and fillcdbID object
235   //
236   TString ocdbPath;
237   Int_t run0=0, run1=0;
238   Int_t version=0, subVersion=0;
239   Bool_t parseStatus = ParseInfoFromOcdbString(ocdbString, ocdbPath, run0,run1,version,subVersion); 
240   if (parseStatus) {
241     AliCDBRunRange runRange(run0,run1);
242     cdbId=AliCDBId(ocdbPath.Data(),runRange,version,subVersion);
243     AliCDBId* id = AliCDBId::MakeFromString(ocdbString);
244     cdbId=*id;
245     delete id;
246   }
247   //
248   return parseStatus;
249 }
250
251 TList  * AliOCDBtoolkit::ConvertListStringToCDBId(const TList *cdbList0){
252   //
253   // Convert input  list of the TObjString to list to AliCDBid 
254   //
255   Int_t entriesList0=cdbList0->GetEntries();
256   TList * array0 = new TList();
257   AliCDBId tmp0;
258   for (Int_t ientry0=0; ientry0<entriesList0; ientry0++){
259     if (cdbList0->At(ientry0)==0) continue;
260     Bool_t isId =  cdbList0->At(ientry0)->IsA()->InheritsFrom("AliCDBId");
261     if (isId){
262       array0->AddLast(cdbList0->At(ientry0));
263     }else{
264       Bool_t isString =  cdbList0->At(ientry0)->IsA()->InheritsFrom("TObjString");
265       if (isString){
266         TObjString* sid0 = dynamic_cast<TObjString*> (cdbList0->At(ientry0));
267         Bool_t status =  ParseInfoFromOcdbString(sid0->String(), tmp0);
268         if (!status) continue;
269         array0->AddLast(new AliCDBId(tmp0));
270       }
271     }
272   }
273   return array0;  
274 }
275
276
277
278 void AliOCDBtoolkit::LoadOCDBFromLog(const char *logName, Int_t verbose){
279   //
280   // Initilaize OCDB
281   // Load OCDB setting as specified in log
282   // Assuming fixed version of the log 
283   // AliCDBManager is initilaized - ocdbMap and ID list can be exported
284   //
285
286   // Parsing/loading sequence:
287   //    0.) SetDefault storage  *** Default Storage URI:
288   //    1.) SetSpecific storage  *** Specific storage
289   //    2.) SetRunNumber  Run number:
290   //    3.) Set used IDs
291   //
292   AliCDBManager * man = AliCDBManager::Instance();
293   //
294   // 0.) SetDefault storage  *** Default Storage URI:
295   // 
296   TString  defaultOCDB = gSystem->GetFromPipe(TString::Format("cat %s| grep \"Storage URI:\"",logName).Data());
297   TObjArray *array = defaultOCDB.Tokenize("\"");
298   man->SetDefaultStorage(array->Last()->GetName());
299   delete array;
300   //
301   // 1.) SetSpecific storage  *** Specific storage
302   //
303   TString  specificStorage  = gSystem->GetFromPipe(TString::Format("cat %s| grep \"Specific storage\"",logName).Data());
304   array = specificStorage.Tokenize("\"");
305   Int_t entries = array->GetEntries();
306   for (Int_t i=1; i<entries-2; i+=4){    
307     // add protection here line shuld be in expected format
308     if ((verbose&2)>0) printf("%s\t%s\n",array->At(i)->GetName(),array->At(i+2)->GetName());    
309     man->SetSpecificStorage(array->At(i)->GetName(),array->At(i+2)->GetName());
310   }
311   delete array;
312   //
313   // 2.) SetRunNumber  Run number:
314   //
315   TString  runLine  = gSystem->GetFromPipe(TString::Format("cat %s| grep \"I-AliCDBManager::Print: Run number =\"",logName).Data());
316   array = runLine.Tokenize("=");
317   Int_t run = 0;
318   if (array->GetEntries()>1) run=atoi(array->At(1)->GetName()); 
319   delete array;
320   man->SetRun(run);  
321   //
322   // 3.) Set used IDs
323   //   
324   TString  ids =   gSystem->GetFromPipe(TString::Format("cat %s| grep I-AliCDB | grep path| grep range | grep version", logName).Data());
325   array= ids.Tokenize("\n");
326   entries = array->GetEntries();
327   //
328   for (Int_t i=0; i<entries; i++){
329     //
330     TString ocdbString = array->At(i)->GetName();
331     TString ocdbEntry;
332     TString ocdbPath;
333     Int_t run0=0, run1=0;
334     Int_t version=0, subVersion=0;
335     Bool_t parseStatus = ParseInfoFromOcdbString(ocdbString, ocdbPath, run0,run1,version,subVersion); 
336     if (!parseStatus) continue;
337     AliCDBRunRange runRange(run0,run1);
338     //
339     if ((verbose&2)!=0) {
340       printf("%s/Run%d_%d_v%d_s%d.root\n",ocdbPath.Data(),run0,run1,version,subVersion); 
341     }
342     try {
343       man->Get(ocdbPath.Data(),runRange,version,subVersion);      
344     } catch(const exception &e){
345       cerr << "OCDB retrieval failed!" << endl;
346       cerr << "Detailes: " << e.what() << endl;
347     }    
348   }  
349   if ((verbose&1)!=0){
350     man->Print();
351     man->GetStorageMap()->Print();
352     man->GetRetrievedIds()->Print(); 
353   }
354 }
355
356 void  AliOCDBtoolkit::SetStorage(const TMap *cdbMap){
357   //
358   // Set storages as speified in the map - TO CHECK.. Should go to the AliCDBmanager if not alreadyhhere
359   //   
360   AliCDBManager * man = AliCDBManager::Instance();  
361   TIter iter(cdbMap->GetTable());
362   TPair* aPair=0;
363   while ((aPair = (TPair*) iter.Next())) {
364     //    aPair->Value();
365     //aPair->Print();
366     if (TString(aPair->GetName())=="default") man->SetDefaultStorage(aPair->Value()->GetName());
367     else
368       man->SetSpecificStorage(aPair->GetName(), aPair->Value()->GetName());
369   }  
370 }
371  
372 void AliOCDBtoolkit::LoadOCDBFromMap(const TMap *cdbMap, const TList *cdbList){
373   //
374   // Initilaize OCDB
375   // Load OCDB setting as specified in maps
376   // Or Do we have already implementation in AliCDBanager?  TO CHECK.. Should go to the AliCDBmanager if not alreadyhhere
377   AliCDBManager * man = AliCDBManager::Instance();  
378   AliOCDBtoolkit::SetStorage(cdbMap);  
379   TIter iter(cdbList);
380   TObjString *ocdbString=0;
381   while (( ocdbString= (TObjString*) iter.Next())) {
382     AliCDBId* cdbId = AliCDBId::MakeFromString(ocdbString->String());
383     try {
384       //      AliCDBEntry * cdbEntry = (AliCDBEntry*) man->Get(*cdbId,kTRUE);
385       man->Get(*cdbId,kTRUE);
386     } catch(const exception &e){
387       cerr << "OCDB retrieval failed!" << endl;
388       cerr << "Detailes: " << e.what() << endl;
389     }   
390   }    
391 }
392
393 void AliOCDBtoolkit::LoadOCDBFromESD(const char *fname){
394   //
395   // Load OCDB setup from the ESD file
396   // 
397   TFile * fesd = TFile::Open(fname);
398   TList *listESD = ((TTree*)fesd->Get("esdTree"))->GetUserInfo();
399   TMap *cdbMapESD= (TMap*)listESD->FindObject("cdbMap");  
400   TList *cdbListESD0= (TList*)listESD->FindObject("cdbList"); // this is list of TObjStrings
401   AliOCDBtoolkit::LoadOCDBFromMap(cdbMapESD, cdbListESD0);
402 }
403
404
405 void AliOCDBtoolkit::MakeDiff(const TMap */*cdbMap0*/, const TList *cdbList0, const TMap */*cdbMap1*/, const TList *cdbList1, Int_t /*verbose*/){
406   //
407   //
408   // Print difference between the 2 ocdb maps
409   // Input:
410   //   maps and list charactireizing OCDB setup  
411   // Output:
412   //   To be decided.
413   //
414   Int_t entriesList0=cdbList0->GetEntries();
415   Int_t entriesList1=cdbList1->GetEntries();
416   //
417   for (Int_t ientry0=0; ientry0<entriesList0; ientry0++){
418     AliCDBId *id0    = dynamic_cast<AliCDBId*> (cdbList0->At(ientry0));
419     AliCDBId *id1=0;
420     for (Int_t ientry1=0; ientry1<entriesList1; ientry1++){
421       AliCDBId *cid1    = dynamic_cast<AliCDBId*> (cdbList1->At(ientry1));
422       //id0.Print();
423       //cid1.Print();
424       if (cid1->GetPath().Contains(id0->GetPath().Data())==0) continue;
425       id1=cid1;
426     }
427     if (!id1) {
428       printf("Missing entry\t");
429       id0->Print();
430       continue;
431     }
432     //   Bool_t isOK=kTRUE;
433     if (id0->GetFirstRun()!= id1->GetFirstRun() ||id0->GetLastRun()!= id1->GetLastRun()){
434       printf("Differrent run range\n");
435       id0->Print();
436       id1->Print();
437     }    
438     if (id0->GetVersion()!= id1->GetVersion() ||id0->GetSubVersion()!= id1->GetSubVersion()){
439       printf("Differrent version\n");
440       id0->Print();
441       id1->Print();
442     }    
443   }
444 }
445
446 void AliOCDBtoolkit::DumpOCDB(const TMap *cdbMap0, const TList *cdbList0, const TString outfile){
447   //
448   // Dump the OCDB configuatation as formated text file 
449   // with following collumns
450   // cdb name  prefix cdb path
451   // OCDB entries are sorted alphabetically
452   // e.g:
453   // TPC/Calib/RecoParam /hera/alice/jwagner/software/aliroot/AliRoot_TPCdev/OCDB/ TPC/Calib/RecoParam/Run0_999999999_v0_s0.root $SIZE_AliCDBEntry_Object $HASH_AliCDBEntry_Object
454   
455   AliCDBManager * man = AliCDBManager::Instance();
456   
457   TList * cdbList = (TList*) cdbList0;   // sorted array
458   cdbList->Sort();
459
460   TIter next(cdbList);
461   AliCDBId *CDBId=0;
462   TString cdbName="";
463   TString cdbPath="";
464   TObjString *ostr;
465   AliCDBEntry *cdbEntry=0;
466   UInt_t hash;
467   TMessage * file;
468   Int_t size; 
469   FILE *ofs = fopen(outfile.Data(),"w");
470   
471   while ((CDBId  =(AliCDBId*) next())){
472     cdbName = CDBId->GetPath();
473     ostr = (TObjString*)cdbMap0->GetValue(cdbName.Data());
474     if(!ostr) ostr = (TObjString*)cdbMap0->GetValue("default");
475     cdbPath = ostr->GetString();
476     if(cdbPath.Contains("local://"))cdbPath=cdbPath(8,cdbPath.Length()).Data();
477     try {
478       cdbEntry = (AliCDBEntry*) man->Get(*CDBId,kTRUE);
479     }catch(const exception &e){
480       cerr << "OCDB retrieval failed!" << endl;
481       cerr << "Detailes: " << e.what() << endl;
482     }  
483     if (!cdbEntry) {
484       printf("Object not avaliable\n");
485       CDBId->Print();
486       continue;
487     }
488     TObject *obj = cdbEntry->GetObject();
489     file = new TMessage(TBuffer::kWrite);
490     file->WriteObject(obj);
491     size = file->Length();
492     if(!obj){
493       fprintf(ofs,"object %s empty!\n",cdbName.Data());
494       continue;
495     }
496     hash = TString::Hash(file->Buffer(),size);
497     fprintf(ofs,"%s\t%s\t%s/Run%d_%d_v%d_s%d.root\t%d\t%u\n",
498            cdbName.Data(),
499            cdbPath.Data(),
500            cdbName.Data(),
501            CDBId->GetFirstRun(),
502            CDBId->GetLastRun(),
503            CDBId->GetVersion(),
504            CDBId->GetSubVersion(),
505            size,
506            hash
507            );
508     //if(!(CDBId->GetPathLevel(0)).Contains("TPC")) continue;
509     //cout << CDBId.ToString() << endl;
510     delete file;
511   }
512   fclose(ofs);
513 }
514
515
516 //====================================================================================================
517 //  Dump object part
518 //==================================================================================================== 
519
520
521
522
523
524 void AliOCDBtoolkit::DumpOCDBFile(const char *finput , const char *foutput, Bool_t dumpMetaData, Bool_t xml){
525   //
526   //  
527   //  DumpOCDBFile("$ALICE_ROOT/OCDB/ITS/Align/Data/Run0_999999999_v0_s0.root", "ITS_Align_Data_Run0_999999999_v0_s0.dump")
528   //
529   if (finput==0) return ;
530   TFile *falignITS  = TFile::Open(finput);
531   AliCDBEntry *entry  = (AliCDBEntry*)falignITS->Get("AliCDBEntry");
532   if (!entry) return; 
533   TObject *obj = ((AliCDBEntry*)falignITS->Get("AliCDBEntry"))->GetObject();  
534
535   //
536   if (!xml){
537     if (dumpMetaData) gROOT->ProcessLine(TString::Format("((TObject*)%p)->Dump(); >%s",entry, foutput).Data());
538     if (!obj) return;
539     gROOT->ProcessLine(TString::Format("AliOCDBtoolkit::DumpObjectRecursive((TObject*)%p); >>%s",obj, foutput).Data());
540   }
541   if (xml){
542     TFile * f = TFile::Open(TString::Format("%s.xml",foutput).Data(),"recreate");
543     if (dumpMetaData) entry->Write("AliCDBEntry");
544     else obj->Write("AliCDBEntry");
545     f->Close();
546   }
547 }
548
549
550
551 void AliOCDBtoolkit::DumpObjectRecursive(TObject *obj){
552   //
553   //
554   //
555   Int_t counterRec=0;
556   printf("==> Dumping object at: %p, name=%s, class=%s)\n", obj, obj->GetName(), (obj->IsA()->GetName()));
557   DumpObjectRecursive(obj, TString(obj->IsA()->GetName())+".",counterRec);
558 }
559  
560 //
561 //
562 //
563 void AliOCDBtoolkit::DumpObjectRecursive(TObject *obj, TString prefix, Int_t &counterRec){
564   //
565   // Recursive dump of the TObject
566   // Dump all basic types and follow pointers to the objects
567   // current limitation:
568   //    a.) clases and structures not derived from TObject not followed (to fix)
569   //    b.) dynamic arrays not followed
570   //    c.) std maps,array ....  not followed
571   //    
572   //
573   if (!obj) return;
574   //
575   // Special case of Collection classes
576   //
577   if (obj->IsA()->InheritsFrom(TCollection::Class())) {
578     TIter myiter((TCollection*)obj);
579     TObject  *arObject=0;
580     Int_t counter=0;
581     while ((arObject = (TObject*)myiter.Next())) {
582       TString prefixArr = TString::Format("%s[%d]",prefix.Data(),counter);
583       DumpObjectRecursive(arObject,prefixArr,counterRec);
584       counter++;
585     } 
586     counterRec++;
587     return;
588   }
589
590   TClass * cl = obj->IsA();
591   if (!(cl->GetListOfRealData())) cl->BuildRealData();
592   TRealData* rd = 0;
593   TIter next(cl->GetListOfRealData());  
594   while ((rd = (TRealData*) next())) {
595     counterRec++;
596     TDataMember* dm = rd->GetDataMember();
597     TDataType* dtype = dm->GetDataType();
598     Int_t offset = rd->GetThisOffset();
599     char* pointer = ((char*) obj) + offset;
600     
601     if (dm->IsaPointer()) {
602       // We have a pointer to an object or a pointer to an array of basic types.
603       TClass* clobj = 0;
604       if (!dm->IsBasic()) {
605         clobj = TClass::GetClass(dm->GetTypeName());
606       }
607       if (clobj) {
608         // We have a pointer to an object.
609         //
610         if (!clobj->InheritsFrom(TObject::Class())) {
611           // It must be a TObject object.
612           continue; 
613         }
614         char** apointer = (char**) pointer;
615         TObject* robj = (TObject*) *apointer;
616         //      
617         if(!robj)
618           printf("M:%s%s\n",prefix.Data(),dm->GetName()); // Missing - 0 pointer
619         else{
620           printf("T:%s\t%s%s\n", clobj->GetName(),prefix.Data(), dm->GetName());
621           TString prefixNew=prefix;
622           prefixNew+=dm->GetName();
623           prefixNew+=".";
624           if (robj!=obj) DumpObjectRecursive(robj,prefixNew,counterRec);  // trivial check 
625           if (robj==obj){
626             printf("R:%s\t%s%s\n",clobj->GetName(),prefix.Data(), dm->GetName());
627           }
628         }
629       }
630     } else if (dm->IsBasic()) {
631       //
632       // Basic data type
633       //
634       const char* index = dm->GetArrayIndex();
635       if (dm->GetArrayDim()==0){
636         printf("B:\t%s%s\t%s\n", prefix.Data(),rd->GetName(), dtype->AsString(pointer));
637       }
638       //
639       // Basic array - fixed length
640       //
641       //      if (dm->GetArrayDim()>0 && strlen(index) != 0){
642       if (dm->GetArrayDim()>0 ){
643         printf("A:\t%s%s\t",prefix.Data(),rd->GetName());
644         Int_t counter=0;
645         for  (Int_t idim=0; idim<dm->GetArrayDim(); idim++){
646           //printf("A:%d\t%d\n", dm->GetArrayDim(),dm->GetMaxIndex(idim));
647           for (Int_t j=0; j<dm->GetMaxIndex(idim); j++){
648             printf("%s\t",dtype->AsString(pointer+dm->GetUnitSize()*counter));
649             counter++;
650             if (counter%5==0) printf("\nA:\t%s%s\t",prefix.Data(),rd->GetName());
651           }
652         }
653         printf("\n");
654       }
655       //
656       // Basic array - dynamic length
657       //
658       if (dm->GetArrayDim()>0 && strlen(index) != 0){
659         //
660         // Dump first only for the moment
661         //  
662         printf("B:\t%s%s\t%s\n",prefix.Data(),rd->GetName(), dtype->AsString(pointer));
663       }
664     } else {
665     }
666   }
667 }  
668
669 //
670 // Small checks to test the TRealData and TDataType
671 //
672
673
674
675 void DumpDataSimple(){
676   //
677   // Dump example for elenatr data types 
678   //
679   TObject *obj = new TVectorD(20);
680   TClass * cl = obj->IsA();
681   if (!cl->GetListOfRealData()) cl->BuildRealData();
682   //
683   TRealData* rd = 0;
684   rd = (TRealData*)(cl->GetListOfRealData()->FindObject("fNrows"));
685   TDataMember* dm = rd->GetDataMember();
686   TDataType* dtype = dm->GetDataType();
687   //
688   Int_t offset = rd->GetThisOffset();
689   char* pointer = ((char*) obj) + offset;
690   printf("%s\n",dtype->AsString(pointer));
691 }
692
693 void DumpDataArray(){
694   //
695   // print array example
696   // 
697   TObject *obj = new TVectorD(20);
698   TClass * cl = obj->IsA();
699   if (!cl->GetListOfRealData()) cl->BuildRealData();
700   TRealData* rd = 0;
701   rd = (TRealData*)(cl->GetListOfRealData()->FindObject("*fElements"));
702   TDataMember* dm = rd->GetDataMember();
703   TDataType* dtype = dm->GetDataType();
704   dtype->Print();
705   //
706   Int_t offset = rd->GetThisOffset();
707   char* pointer = ((char*) obj) + offset; 
708   printf("%s\n",dtype->AsString(pointer));
709 }
710
711 void DumpTObjectArray(){
712   //
713   //
714   //
715   TObjArray *array = new TObjArray(10);
716   for (Int_t i=0; i<10; i++) array->AddLast(new TNamed(Form("n%d",i), Form("n%d",i)));  
717    AliOCDBtoolkit::DumpObjectRecursive(array);
718   //
719   //
720   TObject *obj = array;
721   TClass * cl = obj->IsA();
722   if (!cl->GetListOfRealData()) cl->BuildRealData();
723   TRealData* rd = 0;
724   rd = (TRealData*)(cl->GetListOfRealData()->FindObject("*fCont"));
725   TDataMember* dm = rd->GetDataMember();
726   TDataType* dtype = dm->GetDataType();
727   //
728   Int_t offset = rd->GetThisOffset();
729   char* pointer = ((char*) obj) + offset;
730   char** apointer = (char**) pointer;
731   //we have pointer to pointer here
732   TObject** ppobj = (TObject**) *apointer;
733   (*ppobj)->Print();
734   //
735   TIter myiter(array);
736   TObject  *arObject; 
737   dtype->Print();
738   while ((arObject = (TObject*)myiter.Next())) {
739     AliOCDBtoolkit::DumpObjectRecursive(arObject);
740   } 
741 }
742
743
744 Bool_t AliOCDBtoolkit::AddoptOCDBEntry( const char *finput, const char *output,  Int_t ustartRun, Int_t uendRun){
745   //
746   // Addopt OCDB entry - keeping all of the CDBentry quantities
747   // // Example usage: 
748   //  AliOCDBtoolkit::AddoptOCDBEntry("/cvmfs/alice.gsi.de/alice/simulation/2008/v4-15-Release/Residual/TPC/Calib/ClusterParam/Run127712_130850_v4_s0.root",0,0,AliCDBRunRange::Infinity())
749   TFile * fin = TFile::Open(finput);
750   if (!fin) return kFALSE;
751   AliCDBEntry * entry = (AliCDBEntry*) fin->Get("AliCDBEntry");
752   if (!entry) return kFALSE;
753   
754   AliCDBStorage* pocdbStorage = 0;
755   if (output!=0) AliCDBManager::Instance()->GetStorage(output);
756   else{
757     TString localStorage = "local://"+gSystem->GetFromPipe("pwd")+"/OCDB";
758     pocdbStorage = AliCDBManager::Instance()->GetStorage(localStorage.Data());
759   }
760   //
761   AliCDBId  idIn = entry->GetId();
762   AliCDBMetaData *metaDataIn = entry->GetMetaData();
763
764   AliCDBMetaData *metaData= new AliCDBMetaData();
765   metaData->SetObjectClassName(metaDataIn->GetObjectClassName());
766   metaData->SetResponsible(TString::Format("%s: copy",metaDataIn->GetResponsible()).Data());
767   metaData->SetBeamPeriod(metaDataIn->GetBeamPeriod());
768   //
769   metaData->SetAliRootVersion(metaDataIn->GetAliRootVersion()); //root version
770   metaData->SetComment((TString::Format("%s: copy",metaDataIn->GetComment()).Data()));
771   AliCDBId* id1=NULL;
772   id1=new AliCDBId(idIn.GetPath(), ustartRun, uendRun);
773   pocdbStorage->Put(entry->GetObject(), (*id1), metaData);
774   return kTRUE;
775 }
776
777
778 void AliOCDBtoolkit::MakeSnapshotFromTxt(const TString fInput, const TString outfile, Bool_t singleKeys){
779   //
780   // Make snasphot form the txt file
781   //
782   AliCDBManager * man = AliCDBManager::Instance();
783   LoadOCDBFromList(fInput.Data());
784   man->DumpToSnapshotFile(outfile.Data(), singleKeys);
785
786 }