BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
XtMdcCalib Class Reference

#include <XtMdcCalib.h>

Inheritance diagram for XtMdcCalib:

Public Member Functions

 XtMdcCalib ()
 ~XtMdcCalib ()
void initialize (TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc, IMdcUtilitySvc *mdcUtilitySvc)
void setParam (MdcCalParams &param)
int fillHist (MdcCalEvent *event)
int updateConst (MdcCalibConst *calconst)
void printCut () const
void clear ()
int getHxtKey (int lay, int entr, int lr, int bin) const
Public Member Functions inherited from MdcCalib
 MdcCalib ()
virtual ~MdcCalib ()

Static Public Member Functions

static void fcnXT (Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
static void fcnXtEdge (Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
static double xtFun (double t, double xtpar[])

Static Public Attributes

static std::vector< double > XMEAS
static std::vector< double > TBINCEN
static std::vector< double > ERR
static double Tmax
static double Dmax
static std::vector< double > XMEASED
static std::vector< double > TBINCENED
static std::vector< double > ERRED

Detailed Description

Definition at line 10 of file XtMdcCalib.h.

Constructor & Destructor Documentation

◆ XtMdcCalib()

XtMdcCalib::XtMdcCalib ( )

Definition at line 33 of file XtMdcCalib.cxx.

33 {
34 m_nbin = 50; // m_nbin=50 before 2011-11-16
35 m_nLR = 3;
36 m_nxtpar = 8;
37
38 for ( int lay = 0; lay < 43; lay++ )
39 {
40 if ( lay < 15 ) m_nEntr[lay] = 1;
41 else m_nEntr[lay] = 2; // for entr>0 and entr<0
42 }
43
44 m_tbinw = 10.0;
45 m_fgIni = false;
46}

◆ ~XtMdcCalib()

XtMdcCalib::~XtMdcCalib ( )

Definition at line 48 of file XtMdcCalib.cxx.

48{}

Member Function Documentation

◆ clear()

void XtMdcCalib::clear ( )
virtual

Implements MdcCalib.

Definition at line 50 of file XtMdcCalib.cxx.

50 {
51 cout << "~XtMdcCalib" << endl;
52
53 unsigned int i;
54 for ( i = 0; i < m_hxt.size(); i++ )
55 {
56 delete m_hxt[i];
57 if ( 0 == ( i % 1000 ) ) cout << "delete hxt[" << i << "]" << endl;
58 }
59 delete m_fdXt;
60 m_hxt.clear();
61 m_hxtmap.clear();
62
64}
virtual void clear()=0
Definition MdcCalib.cxx:83

◆ fcnXT()

void XtMdcCalib::fcnXT ( Int_t & npar,
Double_t * gin,
Double_t & f,
Double_t * par,
Int_t iflag )
static

Definition at line 563 of file XtMdcCalib.cxx.

563 {
564 unsigned int i;
565 int ord;
566 Double_t deta;
567 Double_t chisq = 0.0;
568 Double_t dfit;
569
570 for ( i = 0; i < TBINCEN.size(); i++ )
571 {
572 dfit = 0;
573 for ( ord = 0; ord <= 5; ord++ ) { dfit += par[ord] * pow( TBINCEN[i], ord ); }
574 deta = ( dfit - XMEAS[i] ) / ERR[i];
575 chisq += deta * deta;
576 }
577
578 f = chisq;
579}
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
static std::vector< double > ERR
Definition XtMdcCalib.h:31
static std::vector< double > XMEAS
Definition XtMdcCalib.h:29
static std::vector< double > TBINCEN
Definition XtMdcCalib.h:30

Referenced by updateConst().

◆ fcnXtEdge()

void XtMdcCalib::fcnXtEdge ( Int_t & npar,
Double_t * gin,
Double_t & f,
Double_t * par,
Int_t iflag )
static

Definition at line 598 of file XtMdcCalib.cxx.

599 {
600 unsigned int i;
601 Double_t deta;
602 Double_t chisq = 0.0;
603 Double_t dfit;
604
605 for ( i = 0; i < TBINCENED.size(); i++ )
606 {
607 dfit = par[0] * ( TBINCENED[i] - Tmax ) + Dmax;
608 deta = ( dfit - XMEASED[i] ) / ERRED[i];
609 chisq += deta * deta;
610 }
611
612 f = chisq;
613}
static double Dmax
Definition XtMdcCalib.h:34
static double Tmax
Definition XtMdcCalib.h:33
static std::vector< double > ERRED
Definition XtMdcCalib.h:37
static std::vector< double > XMEASED
Definition XtMdcCalib.h:35
static std::vector< double > TBINCENED
Definition XtMdcCalib.h:36

Referenced by updateConst().

◆ fillHist()

int XtMdcCalib::fillHist ( MdcCalEvent * event)
virtual

Implements MdcCalib.

Definition at line 117 of file XtMdcCalib.cxx.

117 {
118 IMessageSvc* msgSvc;
119 Gaudi::svcLocator()->service( "MessageSvc", msgSvc );
120 MsgStream log( msgSvc, "XtMdcCalib" );
121 log << MSG::DEBUG << "XtMdcCalib::fillHist()" << endmsg;
122
123 MdcCalib::fillHist( event );
124
125 // get EsTimeCol
126 bool esCutFg = event->getEsCutFlag();
127 if ( !esCutFg ) return -1;
128
129 int i;
130 int k;
131 int lay;
132 int iLR;
133 int iEntr;
134 int bin;
135 int key;
136 int hid;
137
138 double dr;
139 double dz;
140 double tdr;
141 double doca;
142 double resi;
143 double entrance;
144 int stat;
145
146 double xtpar[MdcCalXtNPars];
147 int nhitlay;
148 bool fgHitLay[MdcCalNLayer];
149
150 if ( !m_fgIni )
151 {
152 for ( lay = 0; lay < MdcCalNLayer; lay++ )
153 {
154 if ( lay < 8 ) m_docaMax[lay] = m_param.maxDocaInner;
155 else m_docaMax[lay] = m_param.maxDocaOuter;
156 }
157
158 for ( lay = 0; lay < MdcCalNLayer; lay++ )
159 {
160 for ( iEntr = 0; iEntr < m_nEntr[lay]; iEntr++ )
161 {
162 for ( iLR = 0; iLR < MdcCalLR; iLR++ )
163 {
164 // m_mdcFunSvc -> getXtpar(lay, iEntr, iLR, xtpar);
165 if ( 0 == iEntr ) m_mdcFunSvc->getXtpar( lay, 8, iLR, xtpar );
166 else if ( 1 == iEntr ) m_mdcFunSvc->getXtpar( lay, 9, iLR, xtpar );
167 m_tm[lay][iEntr][iLR] = xtpar[6];
168 }
169 }
170 }
171 m_fgIni = true;
172 }
173
174 MdcCalRecTrk* rectrk;
175 MdcCalRecHit* rechit;
176
177 int nhit;
178 int ntrk = event->getNTrk();
179 if ( ( ntrk < m_param.nTrkCut[0] ) || ( ntrk > m_param.nTrkCut[1] ) ) return -1;
180 for ( i = 0; i < ntrk; i++ )
181 {
182 rectrk = event->getRecTrk( i );
183 nhit = rectrk->getNHits();
184
185 // dr cut
186 dr = rectrk->getDr();
187 if ( fabs( dr ) > m_param.drCut ) continue;
188
189 // dz cut
190 dz = rectrk->getDz();
191 if ( fabs( dz ) > m_param.dzCut ) continue;
192
193 // momentum cut
194 double p = rectrk->getP();
195 if ( ( fabs( p ) < m_param.pCut[0] ) || ( fabs( p ) > m_param.pCut[1] ) ) continue;
196
197 // select cosmic track with y<0
198 double phi0 = rectrk->getPhi0();
199 double phiTrk = phi0 + CLHEP::halfpi;
200 if ( phiTrk > CLHEP::twopi ) phiTrk -= CLHEP::twopi;
201 if ( m_param.cosmicDwTrk && ( phiTrk < CLHEP::pi ) ) continue;
202
203 for ( lay = 0; lay < MdcCalNLayer; lay++ ) { fgHitLay[lay] = false; }
204
205 for ( k = 0; k < nhit; k++ )
206 {
207 rechit = rectrk->getRecHit( k );
208 lay = rechit->getLayid();
209 fgHitLay[lay] = true;
210 }
211
212 nhitlay = 0;
213 for ( lay = 0; lay < MdcCalNLayer; lay++ )
214 {
215 if ( fgHitLay[lay] ) nhitlay++;
216 }
217 if ( nhitlay < m_param.nHitLayCut ) continue;
218
219 // bool fgNoise = rectrk->getFgNoiseRatio();
220 // if(m_param.noiseCut && (!fgNoise)) continue;
221
222 for ( k = 0; k < nhit; k++ )
223 {
224 rechit = rectrk->getRecHit( k );
225 lay = rechit->getLayid();
226 doca = rechit->getDocaInc();
227 iLR = rechit->getLR();
228 entrance = rechit->getEntra();
229 resi = rechit->getResiInc();
230 tdr = rechit->getTdrift();
231 stat = rechit->getStat();
232
233 if ( 1 != stat ) continue;
234
235 if ( ( fabs( doca ) > m_docaMax[lay] ) || ( fabs( resi ) > m_param.resiCut[lay] ) )
236 { continue; }
237
238 if ( 0 == lay )
239 {
240 if ( !fgHitLay[1] ) continue;
241 }
242 else if ( 42 == lay )
243 {
244 if ( !fgHitLay[41] ) continue;
245 }
246 else
247 {
248 if ( ( !fgHitLay[lay - 1] ) && ( !fgHitLay[lay + 1] ) ) continue;
249 }
250
251 iEntr = m_mdcFunSvc->getXtEntrIndex( entrance );
252 if ( 1 == m_nEntr[lay] ) { iEntr = 0; }
253 else if ( 2 == m_nEntr[lay] )
254 {
255 if ( entrance > 0.0 ) iEntr = 1;
256 else iEntr = 0;
257 }
258
259 bin = (int)( tdr / m_tbinw );
260
261 // left-right separately
262 key = getHxtKey( lay, iEntr, iLR, bin );
263 if ( ( bin < m_nbin ) && ( 1 == m_hxtmap.count( key ) ) )
264 {
265 hid = m_hxtmap[key];
266 m_hxt[hid]->Fill( resi );
267 }
268
269 // left-right in common
270 key = getHxtKey( lay, iEntr, 2, bin );
271 if ( ( bin < m_nbin ) && ( 1 == m_hxtmap.count( key ) ) )
272 {
273 hid = m_hxtmap[key];
274 m_hxt[hid]->Fill( resi );
275 }
276
277 if ( fabs( tdr - m_tm[lay][iEntr][iLR] ) < 5 )
278 {
279 key = getHxtKey( lay, iEntr, iLR, m_nbin );
280 if ( 1 == m_hxtmap.count( key ) )
281 {
282 hid = m_hxtmap[key];
283 m_hxt[hid]->Fill( resi );
284 }
285
286 key = getHxtKey( lay, iEntr, 2, m_nbin );
287 if ( 1 == m_hxtmap.count( key ) )
288 {
289 hid = m_hxtmap[key];
290 m_hxt[hid]->Fill( resi );
291 }
292 }
293
294 } // hit loop
295 } // track loop
296 return 1;
297}
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per bin
Definition FoamA.h:85
const int MdcCalNLayer
Definition MdcCalParams.h:6
const int MdcCalXtNPars
const int MdcCalLR
IMessageSvc * msgSvc()
*************DOUBLE PRECISION m_pi *DOUBLE PRECISION m_HvecTau2 DOUBLE PRECISION m_HvClone2 DOUBLE PRECISION m_gamma1 DOUBLE PRECISION m_gamma2 DOUBLE PRECISION m_thet1 DOUBLE PRECISION m_thet2 INTEGER m_IFPHOT *COMMON c_Taupair $ !Spin Polarimeter vector first Tau $ !Spin Polarimeter vector second Tau $ !Clone Spin Polarimeter vector first Tau $ !Clone Spin Polarimeter vector second Tau $ !Random Euler angle for cloning st tau $ !Random Euler angle for cloning st tau $ !Random Euler angle for cloning st tau $ !Random Euler angle for cloning nd tau $ !Random Euler angle for cloning nd tau $ !Random Euler angle for cloning nd tau $ !phi of HvecTau1 $ !theta of HvecTau1 $ !phi of HvecTau2 $ !theta of HvecTau2 $ !super key
Definition Taupair.h:42
double getDocaInc() const
double getTdrift() const
double getEntra() const
int getLayid() const
double getResiInc() const
int getLR() const
int getStat() const
MdcCalRecHit * getRecHit(int index) const
double getP() const
double getDz() const
int getNHits() const
double getPhi0() const
double getDr() const
virtual int fillHist(MdcCalEvent *event)=0
Definition MdcCalib.cxx:767
int getHxtKey(int lay, int entr, int lr, int bin) const

◆ getHxtKey()

int XtMdcCalib::getHxtKey ( int lay,
int entr,
int lr,
int bin ) const

Definition at line 553 of file XtMdcCalib.cxx.

553 {
554 int key;
555
556 key = ( ( lay << HXTLAYER_INDEX ) & HXTLAYER_MASK ) |
557 ( ( entr << HXTENTRA_INDEX ) & HXTENTRA_MASK ) |
558 ( ( lr << HXTLR_INDEX ) & HXTLR_MASK ) | ( ( bin << HXTBIN_INDEX ) & HXTBIN_MASK );
559
560 return key;
561}

Referenced by fillHist(), initialize(), and updateConst().

◆ initialize()

void XtMdcCalib::initialize ( TObjArray * hlist,
IMdcGeomSvc * mdcGeomSvc,
IMdcCalibFunSvc * mdcFunSvc,
IMdcUtilitySvc * mdcUtilitySvc )
virtual

Implements MdcCalib.

Definition at line 66 of file XtMdcCalib.cxx.

67 {
68 IMessageSvc* msgSvc;
69 Gaudi::svcLocator()->service( "MessageSvc", msgSvc );
70 MsgStream log( msgSvc, "XtMdcCalib" );
71 log << MSG::INFO << "XtMdcCalib::initialize()" << endmsg;
72
73 m_hlist = hlist;
74 m_mdcGeomSvc = mdcGeomSvc;
75 m_mdcFunSvc = mdcFunSvc;
76 m_mdcUtilitySvc = mdcUtilitySvc;
77
78 MdcCalib::initialize( m_hlist, m_mdcGeomSvc, m_mdcFunSvc, m_mdcUtilitySvc );
79
80 int lay;
81 int iLR;
82 int iEntr;
83 int bin;
84 int hid;
85 int key;
86 char hname[200];
87 TH1D* hist;
88
89 m_fdXt = new TFolder( "fdXt", "fdXt" );
90 m_hlist->Add( m_fdXt );
91
92 m_nlayer = m_mdcGeomSvc->getLayerSize();
93
94 hid = 0;
95 for ( lay = 0; lay < m_nlayer; lay++ )
96 {
97 for ( iEntr = 0; iEntr < m_nEntr[lay]; iEntr++ )
98 {
99 for ( iLR = 0; iLR < m_nLR; iLR++ )
100 {
101 for ( bin = 0; bin <= m_nbin; bin++ )
102 {
103 sprintf( hname, "Hxt%02d_E%02d_LR%01d_B%02d", lay, iEntr, iLR, bin );
104 hist = new TH1D( hname, "", 300, -1.5, 1.5 );
105 m_hxt.push_back( hist );
106 m_fdXt->Add( hist );
107
108 key = getHxtKey( lay, iEntr, iLR, bin );
109 m_hxtmap.insert( valType2( key, hid ) );
110 hid++;
111 }
112 }
113 }
114 }
115}
sprintf(cut, "kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_" "pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)
map< int, int >::value_type valType2
virtual void initialize(TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc, IMdcUtilitySvc *mdcUtilitySvc)=0
Definition MdcCalib.cxx:225

◆ printCut()

void XtMdcCalib::printCut ( ) const
virtual

Implements MdcCalib.

Definition at line 299 of file XtMdcCalib.cxx.

virtual void printCut() const =0

◆ setParam()

void XtMdcCalib::setParam ( MdcCalParams & param)
inlinevirtual

Implements MdcCalib.

Definition at line 76 of file XtMdcCalib.h.

76 {
77 MdcCalib::setParam( param );
78 m_param = param;
79}
virtual void setParam(MdcCalParams &param)=0
Definition MdcCalib.h:305

◆ updateConst()

int XtMdcCalib::updateConst ( MdcCalibConst * calconst)
virtual

Implements MdcCalib.

Definition at line 301 of file XtMdcCalib.cxx.

301 {
302 IMessageSvc* msgSvc;
303 Gaudi::svcLocator()->service( "MessageSvc", msgSvc );
304 MsgStream log( msgSvc, "XtMdcCalib" );
305 log << MSG::INFO << "XtMdcCalib::updateConst()" << endmsg;
306
307 MdcCalib::updateConst( calconst );
308
309 int i;
310 int hid;
311 int lay;
312 int iLR;
313 int iEntr;
314 int bin;
315 int ord;
316 int key;
317
318 Stat_t entry;
319 double xtpar;
320 double xterr;
321 double tbcen;
322 double deltx;
323 double xcor;
324 double xerr;
325 double xtini[8];
326 double xtfit[8];
327
328 Int_t ierflg;
329 Int_t istat;
330 Int_t nvpar;
331 Int_t nparx;
332 Double_t fmin;
333 Double_t edm;
334 Double_t errdef;
335 Double_t arglist[10];
336
337 TMinuit* gmxt = new TMinuit( 6 );
338 gmxt->SetPrintLevel( -1 );
339 gmxt->SetFCN( fcnXT );
340 gmxt->SetErrorDef( 1.0 );
341 gmxt->mnparm( 0, "xtpar0", 0, 0.1, 0, 0, ierflg );
342 gmxt->mnparm( 1, "xtpar1", 0, 0.1, 0, 0, ierflg );
343 gmxt->mnparm( 2, "xtpar2", 0, 0.1, 0, 0, ierflg );
344 gmxt->mnparm( 3, "xtpar3", 0, 0.1, 0, 0, ierflg );
345 gmxt->mnparm( 4, "xtpar4", 0, 0.1, 0, 0, ierflg );
346 gmxt->mnparm( 5, "xtpar5", 0, 0.1, 0, 0, ierflg );
347 arglist[0] = 0;
348 gmxt->mnexcm( "SET NOW", arglist, 0, ierflg );
349
350 TMinuit* gmxtEd = new TMinuit( 1 );
351 gmxtEd->SetPrintLevel( -1 );
352 gmxtEd->SetFCN( fcnXtEdge );
353 gmxtEd->SetErrorDef( 1.0 );
354 gmxtEd->mnparm( 0, "xtpar0", 0, 0.1, 0, 0, ierflg );
355 arglist[0] = 0;
356 gmxtEd->mnexcm( "SET NOW", arglist, 0, ierflg );
357
358 ofstream fxtlog( "xtlog" );
359 for ( lay = 0; lay < m_nlayer; lay++ )
360 {
361 if ( 0 == m_param.fgCalib[lay] ) continue;
362
363 for ( iEntr = 0; iEntr < m_nEntr[lay]; iEntr++ )
364 {
365 for ( iLR = 0; iLR < m_nLR; iLR++ )
366 {
367
368 fxtlog << "Layer " << setw( 3 ) << lay << setw( 3 ) << iEntr << setw( 3 ) << iLR
369 << endl;
370
371 for ( ord = 0; ord < m_nxtpar; ord++ )
372 {
373 // xtpar = calconst -> getXtpar(lay, iEntr, iLR, ord);
374 if ( 0 == iEntr ) xtpar = calconst->getXtpar( lay, 8, iLR, ord );
375 else if ( 1 == iEntr ) xtpar = calconst->getXtpar( lay, 9, iLR, ord );
376
377 xtini[ord] = xtpar;
378 xtfit[ord] = xtpar;
379 }
380
381 Tmax = xtini[6];
382
383 for ( bin = 0; bin <= m_nbin; bin++ )
384 {
385 key = getHxtKey( lay, iEntr, iLR, bin );
386 hid = m_hxtmap[key];
387
388 entry = m_hxt[hid]->GetEntries();
389 if ( entry > 100 )
390 {
391 deltx = m_hxt[hid]->GetMean();
392 xerr = m_hxt[hid]->GetRMS();
393 }
394 else { continue; }
395
396 if ( bin < m_nbin ) tbcen = ( (double)bin + 0.5 ) * m_tbinw;
397 else tbcen = m_tm[lay][iEntr][iLR];
398 xcor = xtFun( tbcen, xtini ) - deltx;
399
400 if ( ( tbcen <= Tmax ) || ( bin == m_nbin ) )
401 {
402 TBINCEN.push_back( tbcen );
403 XMEAS.push_back( xcor );
404 ERR.push_back( xerr );
405 }
406 else
407 {
408 TBINCENED.push_back( tbcen );
409 XMEASED.push_back( xcor );
410 ERRED.push_back( xerr );
411 }
412 fxtlog << setw( 3 ) << bin << setw( 15 ) << deltx << setw( 15 ) << xcor << setw( 15 )
413 << tbcen << setw( 15 ) << xerr << endl;
414 } // end of bin loop
415
416 if ( XMEAS.size() < 12 )
417 {
418 TBINCEN.clear();
419 XMEAS.clear();
420 ERR.clear();
421
422 TBINCENED.clear();
423 XMEASED.clear();
424 ERRED.clear();
425
426 continue;
427 }
428
429 for ( ord = 0; ord <= 5; ord++ )
430 {
431 arglist[0] = ord + 1;
432 arglist[1] = xtini[ord];
433 gmxt->mnexcm( "SET PARameter", arglist, 2, ierflg );
434 }
435
436 // fix the xtpar[0] at 0
437 if ( 1 == m_param.fixXtC0 )
438 {
439 arglist[0] = 1;
440 arglist[1] = 0.0;
441 gmxt->mnexcm( "SET PARameter", arglist, 2, ierflg );
442 gmxt->mnexcm( "FIX", arglist, 1, ierflg );
443 }
444
445 arglist[0] = 1000;
446 arglist[1] = 0.1;
447 gmxt->mnexcm( "MIGRAD", arglist, 2, ierflg );
448 gmxt->mnstat( fmin, edm, errdef, nvpar, nparx, istat );
449
450 fxtlog << "Xtpar: " << endl;
451 if ( ( 0 == ierflg ) && ( istat >= 2 ) )
452 {
453 for ( ord = 0; ord <= 5; ord++ )
454 {
455 gmxt->GetParameter( ord, xtpar, xterr );
456 // calconst -> resetXtpar(lay, iEntr, iLR, ord, xtpar);
457 xtfit[ord] = xtpar;
458
459 if ( 1 == m_nEntr[lay] )
460 {
461 for ( i = 0; i < 18; i++ ) calconst->resetXtpar( lay, i, iLR, ord, xtpar );
462 }
463 else if ( 2 == m_nEntr[lay] )
464 {
465 if ( 0 == iEntr )
466 {
467 for ( i = 0; i < 9; i++ ) // entr<0
468 calconst->resetXtpar( lay, i, iLR, ord, xtpar );
469 }
470 else
471 {
472 for ( i = 9; i < 18; i++ ) // entr>0
473 calconst->resetXtpar( lay, i, iLR, ord, xtpar );
474 }
475 }
476 fxtlog << setw( 15 ) << xtpar << setw( 15 ) << xterr << endl;
477 }
478 }
479 else
480 {
481 for ( ord = 0; ord <= 5; ord++ )
482 { fxtlog << setw( 15 ) << xtini[ord] << setw( 15 ) << "0" << endl; }
483 }
484 fxtlog << setw( 15 ) << Tmax << setw( 15 ) << "0" << endl;
485
486 // release the first parameter
487 if ( 1 == m_param.fixXtC0 )
488 {
489 arglist[0] = 1;
490 gmxt->mnexcm( "REL", arglist, 1, ierflg );
491 }
492
493 Dmax = xtFun( Tmax, xtfit );
494
495 if ( XMEASED.size() >= 3 )
496 {
497 // fit xt in the edge area
498 arglist[0] = 1;
499 arglist[1] = xtini[7];
500 gmxtEd->mnexcm( "SET PARameter", arglist, 2, ierflg );
501
502 arglist[0] = 1000;
503 arglist[1] = 0.1;
504 gmxtEd->mnexcm( "MIGRAD", arglist, 2, ierflg );
505 gmxtEd->mnstat( fmin, edm, errdef, nvpar, nparx, istat );
506
507 if ( ( 0 == ierflg ) && ( istat >= 2 ) )
508 {
509 gmxtEd->GetParameter( 0, xtpar, xterr );
510 if ( xtpar < 0.0 ) xtpar = 0.0;
511 if ( m_param.fixXtEdge ) xtpar = xtini[7];
512 // calconst -> resetXtpar(lay, iEntr, iLR, 7, xtpar);
513
514 if ( 1 == m_nEntr[lay] )
515 {
516 for ( i = 0; i < 18; i++ ) calconst->resetXtpar( lay, i, iLR, 7, xtpar );
517 }
518 else if ( 2 == m_nEntr[lay] )
519 {
520 if ( 0 == iEntr )
521 {
522 for ( i = 0; i < 9; i++ ) calconst->resetXtpar( lay, i, iLR, 7, xtpar );
523 }
524 else
525 {
526 for ( i = 9; i < 18; i++ ) calconst->resetXtpar( lay, i, iLR, 7, xtpar );
527 }
528 }
529 fxtlog << setw( 15 ) << xtpar << setw( 15 ) << xterr << endl;
530 }
531 else { fxtlog << setw( 15 ) << xtini[7] << setw( 15 ) << "0" << endl; }
532 }
533 else { fxtlog << setw( 15 ) << xtini[7] << setw( 15 ) << "0" << endl; }
534 fxtlog << "Tm " << setw( 15 ) << Tmax << " Dmax " << setw( 15 ) << Dmax << endl;
535
536 TBINCEN.clear();
537 XMEAS.clear();
538 ERR.clear();
539 TBINCENED.clear();
540 XMEASED.clear();
541 ERRED.clear();
542 } // end of l-r loop
543 } // end of entrance angle loop
544 } // end of layer loop
545
546 fxtlog.close();
547 delete gmxt;
548 delete gmxtEd;
549 cout << "Xt update finished" << endl;
550 return 1;
551}
void resetXtpar(int lay, int entr, int lr, int order, double val)
double getXtpar(int lay, int entr, int lr, int order)
virtual int updateConst(MdcCalibConst *calconst)=0
static double xtFun(double t, double xtpar[])
static void fcnXT(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
static void fcnXtEdge(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)

◆ xtFun()

double XtMdcCalib::xtFun ( double t,
double xtpar[] )
static

Definition at line 581 of file XtMdcCalib.cxx.

581 {
582 int ord;
583 double dist = 0.0;
584 double tmax = xtpar[6];
585
586 if ( t < tmax )
587 {
588 for ( ord = 0; ord <= 5; ord++ ) { dist += xtpar[ord] * pow( t, ord ); }
589 }
590 else
591 {
592 for ( ord = 0; ord <= 5; ord++ ) { dist += xtpar[ord] * pow( tmax, ord ); }
593 dist += xtpar[7] * ( t - tmax );
594 }
595
596 return dist;
597}
int t()
Definition t.c:1

Referenced by updateConst().

Member Data Documentation

◆ Dmax

double XtMdcCalib::Dmax
static

Definition at line 34 of file XtMdcCalib.h.

Referenced by fcnXtEdge(), and updateConst().

◆ ERR

vector< double > XtMdcCalib::ERR
static

Definition at line 31 of file XtMdcCalib.h.

Referenced by fcnXT(), and updateConst().

◆ ERRED

vector< double > XtMdcCalib::ERRED
static

Definition at line 37 of file XtMdcCalib.h.

Referenced by fcnXtEdge(), and updateConst().

◆ TBINCEN

vector< double > XtMdcCalib::TBINCEN
static

Definition at line 30 of file XtMdcCalib.h.

Referenced by fcnXT(), and updateConst().

◆ TBINCENED

vector< double > XtMdcCalib::TBINCENED
static

Definition at line 36 of file XtMdcCalib.h.

Referenced by fcnXtEdge(), and updateConst().

◆ Tmax

double XtMdcCalib::Tmax
static

Definition at line 33 of file XtMdcCalib.h.

Referenced by fcnXtEdge(), and updateConst().

◆ XMEAS

vector< double > XtMdcCalib::XMEAS
static

Definition at line 29 of file XtMdcCalib.h.

Referenced by fcnXT(), and updateConst().

◆ XMEASED

vector< double > XtMdcCalib::XMEASED
static

Definition at line 35 of file XtMdcCalib.h.

Referenced by fcnXtEdge(), and updateConst().


The documentation for this class was generated from the following files: