00001
00002 #include "GaudiKernel/MsgStream.h"
00003 #include "RootMdcCalibDataCnv.h"
00004 #include "CalibData/Mdc/MdcCalibData.h"
00005 #include "CalibData/Mdc/MdcCalStruct.h"
00006 #include "CalibDataSvc/IInstrumentName.h"
00007
00008 #include "TFile.h"
00009 #include "TTree.h"
00010 #include "TDirectory.h"
00011 #include "TObject.h"
00012 #include "TBufferFile.h"
00013 #include "TObjArray.h"
00014
00015 #include "GaudiKernel/CnvFactory.h"
00016 #include "GaudiKernel/IOpaqueAddress.h"
00017 #include "GaudiKernel/DataObject.h"
00018 #include "GaudiKernel/IAddressCreator.h"
00019 #include "GaudiKernel/IDataProviderSvc.h"
00020 #include "GaudiKernel/IConversionSvc.h"
00021 #include "GaudiKernel/GenericAddress.h"
00022
00023 #include "CalibDataSvc/ICalibRootSvc.h"
00024 #include "CalibDataSvc/ICalibMetaCnvSvc.h"
00025
00026
00027 #include "CalibData/CalibModel.h"
00028 using namespace CalibData;
00029
00030
00031
00032
00033
00034 RootMdcCalibDataCnv::RootMdcCalibDataCnv( ISvcLocator* svc) :
00035 RootCalBaseCnv(svc, CLID_Calib_MdcCal) {
00036
00037 }
00038
00039
00040 const CLID& RootMdcCalibDataCnv::objType() const {
00041 return CLID_Calib_MdcCal;
00042 }
00043
00044 const CLID& RootMdcCalibDataCnv::classID() {
00045 return CLID_Calib_MdcCal;
00046 }
00047
00048 StatusCode RootMdcCalibDataCnv::i_createObj(const std::string& fname,
00049 DataObject*& refpObject) {
00050
00051 MsgStream log(msgSvc(), "RootMdcCalibDataCnv");
00052 log<<MSG::DEBUG<<"SetProperty"<<endreq;
00053
00054 StatusCode sc = openRead(fname);
00055 if(!sc)
00056 { log<<MSG::ERROR<<"unable to open files"<<endreq;
00057 }
00058
00059 CalibData::MdcCalibData *tmpObject = new CalibData::MdcCalibData ;
00060
00061 int i;
00062 int nentries;
00063
00064 double xtpar;
00065 int xtkey;
00066 TTree *xttree = (TTree*)m_inFile -> Get("XtTree");
00067 xttree -> SetBranchAddress("xtpar", &xtpar);
00068 xttree -> SetBranchAddress("xtkey", &xtkey);
00069 nentries = xttree -> GetEntries();
00070 for(i=0; i<nentries; i++){
00071 xttree -> GetEntry(i);
00072 tmpObject -> setXtpar(xtkey,xtpar);
00073 }
00074
00075
00076 TObjArray newXtTrees;
00077 if(NULL!=m_inFile->Get("trNewXt00_00_0")){
00078 for(int layid=0; layid<43; layid++){
00079 for(int entr=0; entr<18; entr++){
00080 for(int lr=0; lr<2; lr++){
00081 char newXtTreeName[20];
00082 sprintf(newXtTreeName,"trNewXt%02d_%02d_%d",layid,entr,lr);
00083 TTree* newXtTree = ((TTree*)m_inFile->Get(newXtTreeName));
00084 if(!newXtTree) break;
00085 newXtTrees.Add(newXtTree->CloneTree());
00086 }
00087 }
00088 }
00089 if((43*18*2)==newXtTrees.GetEntries())tmpObject->setNewXtpar(&newXtTrees);
00090 }
00091
00092
00093
00094 TObjArray r2tTrees;
00095 if(NULL!=m_inFile->Get("r2t00")){
00096 for(int layid=0; layid<43; layid++){
00097 char r2tTreeName[20];
00098 sprintf(r2tTreeName,"r2t%02d",layid);
00099 TTree* r2tTree = ((TTree*)m_inFile->Get(r2tTreeName));
00100 if(!r2tTree) break;
00101 r2tTrees.Add(r2tTree->CloneTree());
00102 }
00103 if(43==r2tTrees.GetEntries()) tmpObject->setR2tpar(&r2tTrees);
00104 }
00105
00106
00107 double t0;
00108 double delt0;
00109 TTree *t0tree = (TTree*)m_inFile -> Get("T0Tree");
00110 t0tree -> SetBranchAddress("t0", &t0);
00111 t0tree -> SetBranchAddress("delt0", &delt0);
00112 nentries = t0tree -> GetEntries();
00113 for(i=0; i<nentries; i++){
00114 t0tree -> GetEntry(i);
00115 tmpObject -> setT0(t0);
00116 tmpObject -> setDelT0(delt0);
00117 }
00118
00119
00120 double qtpar0;
00121 double qtpar1;
00122 TTree *qttree = (TTree*)m_inFile -> Get("QtTree");
00123 qttree -> SetBranchAddress("qtpar0", &qtpar0);
00124 qttree -> SetBranchAddress("qtpar1", &qtpar1);
00125 nentries = qttree -> GetEntries();
00126 for(i=0; i<nentries; i++){
00127 qttree -> GetEntry(i);
00128 tmpObject -> setQtpar0(qtpar0);
00129 tmpObject -> setQtpar1(qtpar1);
00130 }
00131
00132
00133 double sdpar;
00134 int sdkey;
00135 TTree *sdtree = (TTree*)m_inFile -> Get("SdTree");
00136 sdtree -> SetBranchAddress("sdpar", &sdpar);
00137 sdtree -> SetBranchAddress("sdkey", &sdkey);
00138 nentries = sdtree -> GetEntries();
00139
00140 for(i=0; i<nentries; i++){
00141 sdtree -> GetEntry(i);
00142 tmpObject -> setSdpar(sdkey,sdpar);
00143 }
00144
00145 refpObject=tmpObject;
00146 return StatusCode::SUCCESS;
00147 }
00148
00149 StatusCode RootMdcCalibDataCnv::createRoot(const std::string& fname,
00150 CalibData::CalibBase1* pTDSObj) {
00151
00152 MsgStream log(msgSvc(), "RootMdcCalibDataCnv");
00153
00154 StatusCode sc = openWrite(fname);
00155 if(!sc)
00156 { log<<MSG::ERROR<<"unable to open files"<<endreq;
00157 }
00158
00159 CalibData::MdcCalibData* tmpObject = dynamic_cast<CalibData::MdcCalibData*>(pTDSObj);
00160 int tmpNo;
00161 int key;
00162 double xtpar;
00163 double t0;
00164 double delt0;
00165 double qtpar[2];
00166 double sdpar;
00167 int i;
00168
00169
00170 TTree *xttree = new TTree("XtTree", "XtTree");
00171 xttree -> Branch("xtkey", &key, "key/I");
00172 xttree -> Branch("xtpar", &xtpar, "xtpar/D");
00173 tmpObject -> setXtBegin();
00174 while( tmpObject -> getNextXtpar(key, xtpar) ){
00175 xttree -> Fill();
00176 }
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199 TTree *t0tree = new TTree("T0Tree", "T0Tree");
00200 t0tree -> Branch("t0", &t0, "t0/D");
00201 t0tree -> Branch("delt0", &delt0, "delt0/D");
00202 tmpNo = tmpObject -> gett0No();
00203 for(i=0; i<tmpNo; i++){
00204 t0 = tmpObject -> getT0(i);
00205 delt0 = tmpObject -> getDelT0(i);
00206 t0tree -> Fill();
00207 }
00208
00209
00210 TTree *qttree = new TTree("QtTree", "QtTree");
00211 qttree -> Branch("qtpar0", &(qtpar[0]), "qtpar0/D");
00212 qttree -> Branch("qtpar1", &(qtpar[1]), "qtpar1/D");
00213 tmpNo = tmpObject -> getqtparNo();
00214 for(i=0; i<tmpNo; i++){
00215 qtpar[0] = tmpObject -> getQtpar0(i);
00216 qtpar[1] = tmpObject -> getQtpar1(i);
00217 qttree -> Fill();
00218 }
00219
00220
00221 TTree *sdtree = new TTree("SdTree", "SdTree");
00222 sdtree -> Branch("sdkey", &key, "key/I");
00223 sdtree -> Branch("sdpar", &sdpar, "sdpar/D");
00224 tmpObject -> setSdBegin();
00225 while( tmpObject -> getNextSdpar(key, sdpar) ){
00226 sdtree -> Fill();
00227 }
00228
00229 xttree -> Write();
00230 t0tree -> Write();
00231 qttree -> Write();
00232 sdtree -> Write();
00233
00234 delete xttree;
00235 delete t0tree;
00236 delete qttree;
00237 delete sdtree;
00238
00239 closeWrite();
00240 log<<MSG::INFO<<"successfully create RootFile"<<endreq;
00241 return sc;
00242 }