#include <Wr2dMdcCalib.h>
Inheritance diagram for Wr2dMdcCalib:
Public Member Functions | |
void | clear () |
void | clear () |
int | fillHist (MdcCalEvent *event) |
int | fillHist (MdcCalEvent *event) |
void | initialize (TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc) |
void | initialize (TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc) |
virtual void | read_file (std::vector< std::string > path) |
virtual void | read_file (std::vector< std::string > path) |
void | setParam (MdcCalParams ¶m) |
void | setParam (MdcCalParams ¶m) |
virtual void | settuple (std::string path) |
virtual void | settuple (std::string path) |
int | updateConst (MdcCalibConst *calconst) |
int | updateConst (MdcCalibConst *calconst) |
Wr2dMdcCalib () | |
Wr2dMdcCalib () | |
~Wr2dMdcCalib () | |
~Wr2dMdcCalib () | |
Static Public Member Functions | |
void | fcnWireParab (Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag) |
void | fcnWireParab (Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag) |
Static Public Attributes | |
bool | fgBIN [MdcCalWrNBin] |
double | xBIN [MdcCalWrNBin] |
double | yBIN [MdcCalWrNBin] |
double | zBIN [MdcCalWrNBin] |
double | zBINERR [MdcCalWrNBin] |
double | zMAX |
double | zMIN |
Private Attributes | |
TFolder * | m_fdWire |
TFolder * | m_fdWire |
TH1F * | m_hl [MdcCalTotCell][MdcCalWrNBin] |
TH1F * | m_hl [MdcCalTotCell][MdcCalWrNBin] |
TObjArray * | m_hlist |
TObjArray * | m_hlist |
TH1F * | m_hr [MdcCalTotCell][MdcCalWrNBin] |
TH1F * | m_hr [MdcCalTotCell][MdcCalWrNBin] |
IMdcCalibFunSvc * | m_mdcFunSvc |
IMdcCalibFunSvc * | m_mdcFunSvc |
IMdcGeomSvc * | m_mdcGeomSvc |
IMdcGeomSvc * | m_mdcGeomSvc |
MdcCalParams | m_param |
double | m_zbinCen [MdcCalNLayer][MdcCalWrNBin] |
double | m_zeast [MdcCalNLayer] |
double | m_zwest [MdcCalNLayer] |
double | m_zwid [MdcCalNLayer] |
|
00027 { 00028 }
|
|
00030 { 00031 }
|
|
|
|
|
|
Implements MdcCalib. |
|
Implements MdcCalib. 00033 { 00034 int bin; 00035 for(int i=0; i<MdcCalTotCell; i++){ 00036 for(bin=0; bin<MdcCalWrNBin; bin++){ 00037 delete m_hl[i][bin]; 00038 delete m_hr[i][bin]; 00039 } 00040 } 00041 delete m_fdWire; 00042 00043 MdcCalib::clear(); 00044 }
|
|
|
|
00362 { 00363 Int_t bin; 00364 Double_t xfit; 00365 Double_t yfit; 00366 00367 double a = 9.47e-06 / (2 * par[4]); 00368 double dx = par[0] - par[2]; 00369 double dy = par[1] - par[3]; 00370 double dz = zMAX - zMIN; 00371 double length = sqrt(dx*dx + dz*dz); 00372 00373 Double_t chisq = 0.0; 00374 Double_t deta; 00375 00376 for(bin=0; bin<MdcCalWrNBin; bin++){ 00377 if( fgBIN[bin] ){ 00378 xfit = (par[0] - par[2]) * (zBIN[bin] - zMIN) / (zMAX - zMIN) + par[2]; 00379 yfit = a*zBIN[bin]*zBIN[bin] + (par[1] - par[3])*zBIN[bin]/length 00380 + 0.5*(par[1] + par[3]) - a*length*length/4.0; 00381 00382 deta = ( (xfit-xBIN[bin])*(xfit-xBIN[bin])+ 00383 (yfit-yBIN[bin])*(yfit-yBIN[bin]) ) / (zBINERR[bin]*zBINERR[bin]); 00384 chisq += deta; 00385 } 00386 } 00387 f = chisq; 00388 }
|
|
Implements MdcCalib. |
|
Implements MdcCalib. 00094 { 00095 IMessageSvc *msgSvc; 00096 Gaudi::svcLocator() -> service("MessageSvc", msgSvc); 00097 MsgStream log(msgSvc, "Wr2dMdcCalib"); 00098 log << MSG::DEBUG << "Wr2dMdcCalib::fillHist()" << endreq; 00099 00100 MdcCalib::fillHist(event); 00101 00102 // get EsTimeCol 00103 bool esCutFg = event->getEsCutFlag(); 00104 if( ! esCutFg ) return -1; 00105 00106 int i; 00107 int k; 00108 int ntrk; 00109 int nhit; 00110 00111 int bin; 00112 int lay; 00113 int cel; 00114 int wir; 00115 int lr; 00116 double dmeas; 00117 double doca; 00118 double residual; 00119 double zhit; 00120 00121 MdcCalRecTrk* rectrk; 00122 MdcCalRecHit* rechit; 00123 00124 ntrk = event->getNTrk(); 00125 log << MSG::DEBUG << "number of tracks: " << ntrk << endreq; 00126 00127 for(i=0; i<ntrk; i++){ 00128 rectrk = event -> getRecTrk(i); 00129 nhit = rectrk -> getNHits(); 00130 00131 log << MSG::DEBUG << "number of hits: " << nhit << endreq; 00132 for(k=0; k<nhit; k++){ 00133 rechit = rectrk -> getRecHit(k); 00134 lay = rechit -> getLayid(); 00135 cel = rechit -> getCellid(); 00136 wir = m_mdcGeomSvc -> Wire(lay, cel) -> Id(); 00137 lr = rechit -> getLR(); 00138 doca = rechit -> getDocaInc(); 00139 dmeas = rechit -> getDmeas(); 00140 residual = rechit -> getResiInc(); 00141 zhit = rechit -> getZhit(); 00142 00143 if((wir<MdcCalTotCell) && (fabs(zhit) < m_zeast[lay])){ 00144 bin = (int)(zhit / m_zwid[lay]); 00145 if( 0 == lr ){ 00146 m_hl[wir][bin] -> Fill(residual); 00147 }else if( 1 == lr ){ 00148 m_hr[wir][bin] -> Fill(residual); 00149 } 00150 }else{ 00151 std::cout << "wir: " << wir << std::endl; 00152 } 00153 } 00154 } 00155 return 1; 00156 }
|
|
Implements MdcCalib. |
|
Implements MdcCalib. 00047 { 00048 IMessageSvc *msgSvc; 00049 Gaudi::svcLocator()->service("MessageSvc", msgSvc); 00050 MsgStream log(msgSvc, "Wr2dMdcCalib"); 00051 log << MSG::INFO << "Wr2dMdcCalib::initialize()" << endreq; 00052 00053 m_hlist = hlist; 00054 m_mdcGeomSvc = mdcGeomSvc; 00055 m_mdcFunSvc = mdcFunSvc; 00056 00057 MdcCalib::initialize(m_hlist, m_mdcGeomSvc, m_mdcFunSvc); 00058 00059 int i; 00060 int bin; 00061 int lay; 00062 int cel; 00063 char hname[200]; 00064 00065 m_fdWire = new TFolder("WireCor", "WireCor"); 00066 m_hlist->Add(m_fdWire); 00067 00068 for(i=0; i<MdcCalTotCell; i++){ 00069 lay = m_mdcGeomSvc -> Wire(i) -> Layer(); 00070 cel = m_mdcGeomSvc -> Wire(i) -> Cell(); 00071 00072 for(bin=0; bin<MdcCalWrNBin; bin++){ 00073 sprintf(hname, "h%04d_L_%02d", i, bin); 00074 m_hl[i][bin] = new TH1F(hname, "", 300, -1.5, 1.5); 00075 m_fdWire->Add(m_hl[i][bin]); 00076 00077 sprintf(hname, "h%04d_R_%02d", i, bin); 00078 m_hr[i][bin] = new TH1F(hname, "", 300, -1.5, 1.5); 00079 m_fdWire->Add(m_hr[i][bin]); 00080 } 00081 } 00082 00083 for(lay=0; lay<MdcCalNLayer; lay++){ 00084 m_zwest[lay] = m_mdcGeomSvc->Wire(lay, 0)->Forward().z(); 00085 m_zeast[lay] = m_mdcGeomSvc->Wire(lay, 0)->Backward().z(); 00086 m_zwid[lay] = (m_zeast[lay] - m_zwest[lay]) / (double)MdcCalWrNBin; 00087 00088 for(bin=0; bin<MdcCalWrNBin; bin++){ 00089 m_zbinCen[lay][bin] = ((double)bin + 0.5) * m_zwid[lay] + m_zwest[lay]; 00090 } 00091 } 00092 }
|
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. 00430 { 00431 MdcCalib::read_file(path); 00432 }
|
|
Implements MdcCalib. |
|
Implements MdcCalib. 00050 { 00051 MdcCalib::setParam(param); 00052 m_param = param; 00053 }
|
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. 00390 { 00391 MdcCalib::settuple(path); 00392 TFile *f1 = new TFile(path.c_str()); 00393 std::cout<<path<<std::endl; 00394 // 00395 // m_mdcGeomSvc = mdcGeomSvc; 00396 // m_mdcFunSvc = mdcFunSvc; 00397 int i; 00398 int bin; 00399 int lay; 00400 int cel; 00401 char hname[200]; 00402 00403 m_fdWire = (TFolder *)f1->Get("WireCor;1"); 00404 00405 for(i=0; i<MdcCalTotCell; i++){ 00406 lay = m_mdcGeomSvc -> Wire(i) -> Layer(); 00407 cel = m_mdcGeomSvc -> Wire(i) -> Cell(); 00408 00409 for(bin=0; bin<MdcCalWrNBin; bin++){ 00410 sprintf(hname, "h%04d_L_%02d", i, bin); 00411 m_hl[i][bin] =(TH1F *)m_fdWire->FindObject(hname); 00412 00413 sprintf(hname, "h%04d_R_%02d", i, bin); 00414 m_hr[i][bin] =(TH1F *)m_fdWire->FindObject(hname); 00415 00416 } 00417 } 00418 00419 for(lay=0; lay<MdcCalNLayer; lay++){ 00420 m_zwest[lay] = m_mdcGeomSvc->Wire(lay, 0)->Forward().z(); 00421 m_zeast[lay] = m_mdcGeomSvc->Wire(lay, 0)->Backward().z(); 00422 m_zwid[lay] = (m_zeast[lay] - m_zwest[lay]) / (double)MdcCalWrNBin; 00423 00424 for(bin=0; bin<MdcCalWrNBin; bin++){ 00425 m_zbinCen[lay][bin] = ((double)bin + 0.5) * m_zwid[lay] + m_zwest[lay]; 00426 } 00427 } 00428 }
|
|
Implements MdcCalib. |
|
Implements MdcCalib. 00158 { 00159 IMessageSvc *msgSvc; 00160 Gaudi::svcLocator() -> service("MessageSvc", msgSvc); 00161 MsgStream log(msgSvc, "Wr2dMdcCalib"); 00162 log << MSG::INFO << "Wr2dMdcCalib::updateConst()" << endreq; 00163 00164 int i; 00165 int k; 00166 int bin; 00167 int lay; 00168 int cel; 00169 double dw; 00170 00171 // minuit for wires 00172 Int_t ierflg; 00173 Int_t istat; 00174 Int_t nvpar; 00175 Int_t nparx; 00176 Double_t fmin; 00177 Double_t edm; 00178 Double_t errdef; 00179 Double_t arglist[10]; 00180 00181 TMinuit *gmtrw = new TMinuit(5); 00182 gmtrw -> SetPrintLevel(-1); 00183 gmtrw -> SetFCN(fcnWireParab); 00184 gmtrw -> SetErrorDef(1.0); 00185 gmtrw -> mnparm(0, "xf", 0.0, 0.01, 0, 0, ierflg); 00186 gmtrw -> mnparm(1, "yf", 0.0, 0.01, 0, 0, ierflg); 00187 gmtrw -> mnparm(2, "xb", 0.0, 0.01, 0, 0, ierflg); 00188 gmtrw -> mnparm(3, "yb", 0.0, 0.01, 0, 0, ierflg); 00189 gmtrw -> mnparm(4, "ten", 0.0, 0.1, 0, 0, ierflg); 00190 arglist[0] = 0; 00191 gmtrw -> mnexcm("Set NOW", arglist, 0, ierflg); 00192 00193 double a; 00194 double dx; 00195 double dy; 00196 double dz; 00197 double length; 00198 double ten[] = {15, 15, 15, 16, 16, 17, 17, 18, 14, 14, 00199 19, 19, 24, 24, 31, 31, 37, 37, 45, 45, 00200 46, 47, 47, 47, 47, 48, 48, 48, 48, 49, 00201 49, 49, 49, 50, 50, 50, 51, 51, 51, 52, 00202 52, 52, 52}; 00203 double wpar[5]; 00204 double wparErr[5]; 00205 double sinphiw; 00206 double cosphiw; 00207 double deldw; 00208 double delxw; 00209 double delyw; 00210 00211 int nFit; 00212 Stat_t entryL; 00213 Stat_t entryR; 00214 Double_t hcen; 00215 Double_t cenL[MdcCalWrNBin]; 00216 Double_t errL[MdcCalWrNBin]; 00217 Double_t cenR[MdcCalWrNBin]; 00218 Double_t errR[MdcCalWrNBin]; 00219 00220 ofstream fwlog("logWireCor.txt"); 00221 ofstream fwpc("wireposCor.txt"); 00222 ofstream fwpcErr("wireposErr.txt"); 00223 00224 fwpc << setw(5) << "wireId" << setw(15) << "dx_East(mm)" 00225 << setw(15) << "dy_East(mm)" << setw(15) << "dz_East(mm)" 00226 << setw(15) << "dx_West(mm)" << setw(15) << "dy_West(mm)" 00227 << setw(15) << "dz_West(mm)" << endl; 00228 00229 for(i=0; i<MdcCalTotCell; i++){ 00230 for(k=0; k<5; k++) wparErr[k] = 999.0; 00231 nFit = 0; 00232 for(bin=0; bin<MdcCalWrNBin; bin++){ 00233 entryL = m_hl[i][bin] -> GetEntries(); 00234 entryR = m_hr[i][bin] -> GetEntries(); 00235 if((entryL < 100) && (entryR < 100)){ 00236 fgBIN[bin] = false; 00237 continue; 00238 } else{ 00239 fgBIN[bin] = true; 00240 } 00241 nFit++; 00242 00243 if(1 == m_param.fgCalib[lay]){ 00244 hcen = m_hl[i][bin] -> GetMean(); 00245 if(entryL > 500){ 00246 m_hl[i][bin] -> Fit("gaus", "Q", "", hcen-1.5, hcen+1.5); 00247 cenL[bin] = m_hl[i][bin]->GetFunction("gaus")->GetParameter(1); 00248 errL[bin] = m_hl[i][bin]->GetFunction("gaus")->GetParameter(2); 00249 } 00250 else{ 00251 cenL[bin] = hcen; 00252 errL[bin] = m_hl[i][bin] -> GetRMS(); 00253 } 00254 00255 hcen = m_hr[i][bin] -> GetMean(); 00256 if(entryR > 500){ 00257 m_hr[i][bin] -> Fit("gaus", "Q", "", hcen-1.5, hcen+1.5); 00258 cenR[bin] = m_hr[i][bin]->GetFunction("gaus")->GetParameter(1); 00259 errR[bin] = m_hr[i][bin]->GetFunction("gaus")->GetParameter(2); 00260 } 00261 else{ 00262 cenR[bin] = hcen; 00263 errR[bin] = m_hr[i][bin] -> GetRMS(); 00264 } 00265 } else{ 00266 cenL[bin] = 0.0; 00267 errL[bin] = 0.15; 00268 cenR[bin] = 0.0; 00269 errR[bin] = 0.15; 00270 } 00271 } 00272 if(nFit < 8) continue; 00273 00274 lay = m_mdcGeomSvc -> Wire(i) -> Layer(); 00275 cel = m_mdcGeomSvc -> Wire(i) -> Cell(); 00276 zMIN = m_zwest[lay]; 00277 zMAX = m_zeast[lay]; 00278 wpar[0] = m_mdcGeomSvc->Wire(lay, 0)->Backward().x(); 00279 wpar[1] = m_mdcGeomSvc->Wire(lay, 0)->Backward().y(); 00280 wpar[2] = m_mdcGeomSvc->Wire(lay, 0)->Forward().x(); 00281 wpar[3] = m_mdcGeomSvc->Wire(lay, 0)->Forward().y(); 00282 wpar[4] = ten[lay]; 00283 00284 a = 9.47e-06 / (2 * wpar[4]); 00285 dx = wpar[0] - wpar[2]; 00286 dy = wpar[1] - wpar[3]; 00287 dz = zMAX - zMIN; 00288 length = sqrt(dx*dx + dz*dz); 00289 00290 for(bin=0; bin<MdcCalWrNBin; bin++){ 00291 zBIN[bin] = m_zbinCen[lay][bin]; 00292 xBIN[bin] = (wpar[0]-wpar[2])*(zBIN[bin]-zMIN)/(zMAX-zMIN) + wpar[2]; 00293 yBIN[bin] = a*zBIN[bin]*zBIN[bin] + (wpar[1]-wpar[3])*zBIN[bin]/length 00294 + 0.5*(wpar[1]+wpar[3]) - a*length*length/4.0; 00295 00296 sinphiw = yBIN[bin] / sqrt(xBIN[bin]*xBIN[bin] + yBIN[bin]*yBIN[bin]); 00297 cosphiw = xBIN[bin] / sqrt(xBIN[bin]*xBIN[bin] + yBIN[bin]*yBIN[bin]); 00298 00299 deldw = - (cenL[bin]-cenR[bin])/2.0; 00300 delxw = -deldw * sinphiw; 00301 delyw = deldw * cosphiw; 00302 00303 fwlog << setw(3) << lay << setw(4) << cel << setw(5) << i 00304 << setw(4) << bin << setw(15) << zBIN[bin] 00305 << setw(15) << deldw << setw(15) << delxw 00306 << setw(15) << delyw << endl; 00307 00308 xBIN[bin] += delxw; 00309 yBIN[bin] += delyw; 00310 00311 zBINERR[bin] = sqrt((errL[bin]*errL[bin]) + (errR[bin]*errR[bin]))/2.0; 00312 } 00313 00314 arglist[0] = 1; 00315 arglist[1] = wpar[0]; 00316 gmtrw -> mnexcm("SET PARameter", arglist, 2, ierflg); 00317 arglist[0] = 2; 00318 arglist[1] = wpar[1]; 00319 gmtrw -> mnexcm("SET PARameter", arglist, 2, ierflg); 00320 arglist[0] = 3; 00321 arglist[1] = wpar[2]; 00322 gmtrw -> mnexcm("SET PARameter", arglist, 2, ierflg); 00323 arglist[0] = 4; 00324 arglist[1] = wpar[3]; 00325 gmtrw -> mnexcm("SET PARameter", arglist, 2, ierflg); 00326 arglist[0] = 5; 00327 arglist[1] = wpar[4]; 00328 gmtrw -> mnexcm("SET PARameter", arglist, 2, ierflg); 00329 gmtrw -> mnexcm("FIX", arglist, 1, ierflg); 00330 00331 arglist[0] = 2000; 00332 arglist[1] = 0.1; 00333 gmtrw -> mnexcm("MIGRAD", arglist, 2, ierflg); 00334 gmtrw -> mnstat(fmin, edm, errdef, nvpar, nparx, istat); 00335 00336 if( (0==ierflg) && (3==istat) ){ 00337 gmtrw -> GetParameter(0, wpar[0], wparErr[0]); 00338 gmtrw -> GetParameter(1, wpar[1], wparErr[1]); 00339 gmtrw -> GetParameter(2, wpar[2], wparErr[2]); 00340 gmtrw -> GetParameter(3, wpar[3], wparErr[3]); 00341 gmtrw -> GetParameter(4, wpar[4], wparErr[4]); 00342 } 00343 gmtrw -> mnexcm("RELease", arglist, 1, ierflg); 00344 00345 fwlog << setw(5) << i << setw(15) << wpar[0] << setw(15) << wpar[1] 00346 << setw(15) << wpar[2] << setw(15) << wpar[3] << endl; 00347 fwpc << setw(5) << i << setw(15) << wpar[0] << setw(15) << wpar[1] 00348 << setw(15) << "0" << setw(15) << wpar[2] 00349 << setw(15) << wpar[3] << setw(15) << "0" << endl; 00350 fwpcErr << setw(5) << i << setw(15) << wparErr[0] << setw(15) << wparErr[1] 00351 << setw(15) << wparErr[2] << setw(15) << wparErr[3] << endl; 00352 } 00353 fwlog.close(); 00354 fwpc.close(); 00355 fwpcErr.close(); 00356 00357 delete gmtrw; 00358 return 1; 00359 }
|
|
|
|
|
|
|
|
|
|
|
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. |
|
|
|
|
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. |
|
Reimplemented from MdcCalib. |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|