00001 #include <iostream>
00002 #include <fstream>
00003 #include <iomanip>
00004 #include <string>
00005 #include <cstring>
00006 #include <vector>
00007
00008 #include "TROOT.h"
00009 #include "TFile.h"
00010 #include "TTree.h"
00011 #include "TFolder.h"
00012 #include "TProfile.h"
00013 #include "TObjArray.h"
00014 #include "TList.h"
00015 #include "TSpline.h"
00016 #include "TPostScript.h"
00017 #include "TLatex.h"
00018 #include "TCanvas.h"
00019 #include "TStyle.h"
00020
00021 #include "include/MdcCalibConst.h"
00022 #include "include/MdcCosGeom.h"
00023 #include "include/fun.h"
00024 #include "include/IniCalib.h"
00025 #include "include/PreXtCalib.h"
00026 #include "include/PreT0Calib.h"
00027 #include "include/XtCalib.h"
00028 #include "include/GrXtCalib.h"
00029 #include "include/XtInteCalib.h"
00030 #include "include/T0Calib.h"
00031 #include "include/QtCalib.h"
00032
00033 using namespace std;
00034
00035 int main(int argc, char* argv[]){
00036 char* jobname;
00037 char* strcal;
00038 int fgCal = 1;
00039 if(argc>2){
00040 jobname = argv[1];
00041 strcal = argv[2];
00042 sscanf(strcal, "%d", &fgCal);
00043 } else if(argc>1){
00044 jobname = argv[1];
00045 }else{
00046 cout << "bad argument" << endl;
00047 return -1;
00048 }
00049 if(fgCal <= 0) cout << "do not calibrate " << endl;
00050
00051 string path = "";
00052 string strJob = jobname;
00053 cout << "strJob: " << strJob << endl;
00054 string::size_type ilast = strJob.find_last_of("/");
00055 if(string::npos != ilast){
00056 path = strJob.substr(0, ilast);
00057 }
00058
00059 int calType;
00060 string constname;
00061 string confname;
00062 string str;
00063 string strtmp;
00064 ifstream fjob(jobname);
00065 if( ! fjob.is_open() ){
00066 cout << "ERROR: can not read jobOption: " << jobname << endl;
00067 return 0;
00068 } else{
00069 cout << "Open jobOption: " << jobname << endl;
00070 while( getline(fjob, str) ){
00071 if(str.find("//", 0) != string::npos){
00072 continue;
00073 } else if( str.find("CalibRootCnvSvc.Mdcrootfile", 0) != string::npos ){
00074 string::size_type i1 = str.find_first_of("\"");
00075 string::size_type i2 = str.find_last_of("\"");
00076 constname = str.substr(i1+1, i2-i1-1);
00077 } else if(str.find("MdcCalibAlg.ConfigFile", 0) != string::npos){
00078 string::size_type i1 = str.find_first_of("\"");
00079 string::size_type i2 = str.find_last_of("\"");
00080 confname = str.substr(i1+1, i2-i1-1);
00081 } else if(str.find("MdcCalibAlg.MdcCalFlg", 0) != string::npos){
00082 string::size_type i1 = str.find_first_of("=");
00083 string::size_type i2 = str.find_last_of(";");
00084 strtmp = str.substr(i1+1, i2-i1-1);
00085 sscanf(strtmp.c_str(), "%d", &calType);
00086 }
00087 }
00088 }
00089
00090 MdcCosGeom* pGeom = 0;
00091 pGeom = new MdcCosGeom("/home/bes/wulh/document/wireconf.txt", "/home/bes/wulh/calibConst/MdcAlignPar_ini.txt" );
00092 pGeom -> initialize(0.0);
00093
00094 int lay;
00095 for(lay=0; lay<15; lay++) gNEntr[lay] = 1;
00096 for(lay=15; lay<43; lay++) gNEntr[lay] = 2;
00097
00098
00099 string strconfig;
00100 string strcomment;
00101 ifstream fconfig(confname.c_str());
00102 if( ! fconfig.is_open() ){
00103 cout << "ERROR: can not read config file" << endl;
00104 return 0;
00105 } else{
00106 cout << "Open config file: " << confname << endl;
00107 while( fconfig >> strconfig ){
00108 if('#' == strconfig[0]){
00109 getline(fconfig, strcomment);
00110 } else if("TimeShift" == strconfig){
00111 fconfig >> gTimeShift;
00112 } else if("TesMin" == strconfig){
00113 fconfig >> gTesMin;
00114 } else if("TesMax" == strconfig){
00115 fconfig >> gTesMax;
00116 } else if("FlagIniCalConst" == strconfig){
00117 fconfig >> gFgIniCalConst;
00118 } else if("FlagUpdateTmInPreT0" == strconfig){
00119 fconfig >> gPreT0SetTm;
00120 } else if("InitT0" == strconfig){
00121 fconfig >> gInitT0;
00122 } else if("T0Shift" == strconfig){
00123 fconfig >> gT0Shift;
00124 } else if("TminFitChindf" == strconfig){
00125 fconfig >> gTminFitChindf;
00126 } else if("TmaxFitChindf" == strconfig){
00127 fconfig >> gTmaxFitChindf;
00128 } else if("ResidualType" == strconfig){
00129 fconfig >> gResiType;
00130 } else if("UpdateSigma" == strconfig){
00131 fconfig >> gCalSigma;
00132 } else if("FixXtC0" == strconfig){
00133 fconfig >> gFixXtC0;
00134 } else if("RawTimeFitRange" == strconfig){
00135 for(lay=0; lay<NLAYER; lay++){
00136 fconfig >> strtmp
00137 >> gFgCalib[lay]
00138 >> gTminFitRange[lay][0]
00139 >> gTminFitRange[lay][1]
00140 >> gTmaxFitRange[lay][0]
00141 >> gTmaxFitRange[lay][1]
00142 >> gInitTm[lay]
00143 >> strtmp
00144 >> gQmin[lay]
00145 >> gQmax[lay];
00146 }
00147 }
00148 }
00149 }
00150 fconfig.close();
00151
00152 TObjArray* hlist = new TObjArray(0);
00153 CalibBase* pcal;
00154 if(0 == calType) pcal = new IniCalib();
00155 else if(1 == calType) pcal = new PreXtCalib();
00156 else if(2 == calType) pcal = new PreT0Calib();
00157 else if(3 == calType) pcal = new XtCalib();
00158 else if(4 == calType) pcal = new GrXtCalib();
00159 else if(9 == calType) pcal = new XtInteCalib();
00160 else if(5 == calType) pcal = new T0Calib();
00161 else if(8 == calType) pcal = new QtCalib();
00162 else {cout << "Error CalibType" << endl; return 0;}
00163 pcal->init(hlist, pGeom);
00164
00165 vector<string> fhistname = getHistList();
00166 if(0==fhistname.size()){
00167 cout << "hist file path: " << path << endl;
00168 fhistname = getHistList(path);
00169 }
00170 for(unsigned nf=0; nf<fhistname.size(); nf++){
00171 TFile* fin = new TFile(fhistname[nf].c_str());
00172 if(!fin->IsOpen()){
00173 continue;
00174 } else{
00175 cout << "merge hist file " << nf << ": " << fhistname[nf] << endl;
00176 pcal->mergeHist(fin);
00177 fin->Close();
00178 }
00179 }
00180
00181
00182 MdcCalibConst* calconst = new MdcCalibConst();
00183
00184 if(constname==""){constname = "/home/bes/wulh/calibConst/MdcCalibConst_11397_psi3770_652b_v1.root";}
00185 TFile fconst(constname.c_str());
00186 if(!fconst.IsOpen()){
00187 cout << "ERROR: " << constname << " does not exist" << endl;
00188 return 0;
00189 }
00190 cout << "Open calib const file " << constname << endl;
00191 int key;
00192 double xtpar;
00193 int entry;
00194 TTree* xttree = (TTree*)fconst.Get("XtTree");
00195 xttree -> SetBranchAddress("xtkey", &key);
00196 xttree -> SetBranchAddress("xtpar", &xtpar);
00197 entry = (int)xttree -> GetEntries();
00198 for(int i=0; i<entry; i++){
00199 xttree -> GetEntry(i);
00200 calconst->fillXtpar(key, xtpar);
00201 }
00202
00203 double t0;
00204 double delt0;
00205 TTree* t0tree = (TTree*)fconst.Get("T0Tree");
00206 t0tree -> SetBranchAddress("t0", &t0);
00207 t0tree -> SetBranchAddress("delt0", &delt0);
00208 entry = (int)t0tree -> GetEntries();
00209 for(int i=0; i<entry; i++){
00210 t0tree -> GetEntry(i);
00211 calconst->fillT0(t0);
00212 calconst->fillDelT0(delt0);
00213 }
00214
00215 double qtpar0;
00216 double qtpar1;
00217 TTree* qttree = (TTree*)fconst.Get("QtTree");
00218 qttree -> SetBranchAddress("qtpar0", &qtpar0);
00219 qttree -> SetBranchAddress("qtpar1", &qtpar1);
00220 entry = (int)qttree -> GetEntries();
00221 for(int i=0; i<entry; i++){
00222 qttree -> GetEntry(i);
00223 calconst->fillQtpar0(qtpar0);
00224 calconst->fillQtpar1(qtpar1);
00225 }
00226
00227 double sdpar;
00228 TTree* sdtree = (TTree*)fconst.Get("SdTree");
00229 sdtree -> SetBranchAddress("sdkey", &key);
00230 sdtree -> SetBranchAddress("sdpar", &sdpar);
00231 entry = sdtree -> GetEntries();
00232 for(int i=0; i<entry; i++){
00233 sdtree -> GetEntry(i);
00234 calconst->fillSdpar(key, sdpar);
00235 }
00236
00237 TObjArray* newXtList = new TObjArray(0);
00238 TList* list = fconst.GetListOfKeys();
00239 int listSize = (int)(list->GetSize());
00240 cout << "Number of trees in old calib const file: " << listSize << endl;
00241 if(listSize > 4){
00242 gROOT->cd();
00243 for (int i = 0; i<listSize; i++) {
00244 TTree* tree = (TTree*)fconst.Get(list->At(i)->GetName());
00245 string strName(tree->GetName());
00246 if(string::npos != strName.find("trNewXt")){
00247 TTree* t2 = tree->CopyTree("");
00248 newXtList->Add(t2);
00249 }
00250 }
00251 }
00252 fconst.Close();
00253
00254 gStyle -> SetCanvasBorderMode(0);
00255 gStyle -> SetCanvasColor(10);
00256 gStyle -> SetOptFit(0011);
00257 gStyle -> SetStatColor(10);
00258 gStyle -> SetStatBorderSize(1);
00259 gStyle -> SetStatFont(42);
00260 gStyle -> SetStatFontSize(0.04);
00261 gStyle -> SetStatX(0.9);
00262 gStyle -> SetStatY(0.9);
00263 gStyle -> SetPadColor(10);
00264 gStyle -> SetFuncColor(2);
00265
00266 TObjArray* r2tList = new TObjArray(0);
00267
00268
00269 pcal->calib(calconst, newXtList, r2tList);
00270
00271 TFile fhist("histall.root", "recreate");
00272 fhist.cd();
00273 hlist->Write();
00274 fhist.Close();
00275
00276
00277 if(fgCal > 0) writeConst(calconst, newXtList, r2tList);
00278
00279 return 0;
00280 }
00281
