source: trunk/Cosy/main/MBending.cc@ 19372

Last change on this file since 19372 was 8376, checked in by tbretz, 19 years ago
*** empty log message ***
File size: 18.4 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Thomas Bretz, 2003 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2000-2003
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MBending
28//
29// Double_t fIe ; // [rad] Index Error in Elevation
30// Double_t fIa ; // [rad] Index Error in Azimuth
31//
32// Double_t fFlop ; // [rad] Vertical Sag
33// * do not use if not data: Zd<0
34//
35// Double_t fNpae ; // [rad] Az-El Nonperpendicularity
36//
37// Double_t fCa ; // [rad] Left-Right Collimation Error
38//
39// Double_t fAn ; // [rad] Azimuth Axis Misalignment (N-S)
40// Double_t fAw ; // [rad] Azimuth Axis Misalignment (E-W)
41//
42// Double_t fTf ; // [rad] Tube fluxture (sin)
43// * same as ecec if no data: Zd<0
44// Double_t fTx ; // [rad] Tube fluxture (tan)
45// * do not use with NPAE if no data: Zd<0
46//
47// Double_t fNrx ; // [rad] Nasmyth rotator displacement, horizontan
48// Double_t fNry ; // [rad] Nasmyth rotator displacement, vertical
49//
50// Double_t fCrx ; // [rad] Alt/Az Coude Displacement (N-S)
51// Double_t fCry ; // [rad] Alt/Az Coude Displacement (E-W)
52//
53// Double_t fEces ; // [rad] Elevation Centering Error (sin)
54// Double_t fAces ; // [rad] Azimuth Centering Error (sin)
55// Double_t fEcec ; // [rad] Elevation Centering Error (cos)
56// Double_t fAcec ; // [rad] Azimuth Centering Error (cos)
57//
58////////////////////////////////////////////////////////////////////////////
59#include "MBending.h"
60
61#include <fstream>
62#include <iomanip>
63#include <TVector3.h>
64
65#include <TMinuit.h>
66
67#include "MTime.h"
68
69ClassImp(MBending);
70
71using namespace std;
72
73#undef DEBUG
74//#define DEBUG(txt) txt
75#define DEBUG(txt)
76
77const Int_t MBending::fNumPar=19;
78
79void MBending::Init()
80{
81 fCoeff = new Double_t*[fNumPar];
82 fName = new TString[fNumPar];
83 fDescr = new TString[fNumPar];
84
85 fCoeff[ 0] = &fIa; fName[ 0] = "IA";
86 fCoeff[ 1] = &fIe; fName[ 1] = "IE";
87 fCoeff[ 2] = &fFlop; fName[ 2] = "FLOP";
88 fCoeff[ 3] = &fAn; fName[ 3] = "AN";
89 fCoeff[ 4] = &fAw; fName[ 4] = "AW";
90 fCoeff[ 5] = &fNpae; fName[ 5] = "NPAE";
91 fCoeff[ 6] = &fCa; fName[ 6] = "CA";
92 fCoeff[ 7] = &fTf; fName[ 7] = "TF";
93 fCoeff[ 8] = &fTx; fName[ 8] = "TX";
94 fCoeff[ 9] = &fEces; fName[ 9] = "ECES";
95 fCoeff[10] = &fAces; fName[10] = "ACES";
96 fCoeff[11] = &fEcec; fName[11] = "ECEC";
97 fCoeff[12] = &fAcec; fName[12] = "ACEC";
98 fCoeff[13] = &fNrx; fName[13] = "NRX";
99 fCoeff[14] = &fNry; fName[14] = "NRY";
100 fCoeff[15] = &fCrx; fName[15] = "CRX";
101 fCoeff[16] = &fCry; fName[16] = "CRY";
102 fCoeff[17] = &fMagic1; fName[17] = "MAGIC1";
103 fCoeff[18] = &fMagic2; fName[18] = "MAGIC2";
104
105 fDescr[ 0] = "Index Error Azimuth";
106 fDescr[ 1] = "Index Error Zenith Distance";
107 fDescr[ 2] = "Vertical Sag";
108 fDescr[ 3] = "Azimuth Axis Misalignment (N-S)";
109 fDescr[ 4] = "Azimuth Axis Misalignment (E-W)";
110 fDescr[ 5] = "Az-El Nonperpendicularity";
111 fDescr[ 6] = "Left-Right Collimation Error";
112 fDescr[ 7] = "Tube fluxture (sin)";
113 fDescr[ 8] = "Tube fluxture (tan)";
114 fDescr[ 9] = "Elevation Centering Error (sin)";
115 fDescr[10] = "Azimuth Centering Error (sin)";
116 fDescr[11] = "Elevation Centering Error (cos)";
117 fDescr[12] = "Azimuth Centering Error (cos)";
118 fDescr[13] = "Nasmyth rotator displacement (horizontal)";
119 fDescr[14] = "Nasmyth rotator displacement (vertical)";
120 fDescr[15] = "Alt/Az Coude Displacement (N-S)";
121 fDescr[16] = "Alt/Az Coude Displacement (E-W)";
122 fDescr[17] = "n/a <ZA Hysteresis>";
123 fDescr[18] = "n/a";
124}
125void MBending::Reset()
126{
127 Clear();
128}
129
130void MBending::Load(const char *name)
131{
132 /*
133 ! MMT 1987 July 8
134 ! T 36 7.3622 41.448 -0.0481
135 ! IA -37.5465 20.80602
136 ! IE -13.9180 1.25217
137 ! NPAE +7.0751 26.44763
138 ! CA -6.9149 32.05358
139 ! AN +0.5053 1.40956
140 ! AW -2.2016 1.37480
141 ! END
142 */
143
144 ifstream fin(name);
145 if (!fin)
146 {
147 cout << "Error: Cannot open file '" << name << "'" << endl;
148 return;
149 }
150
151 char c;
152 while (fin && fin.get()!='\n');
153 fin >> c;
154
155 if (c!='S' && c!='s')
156 {
157 cout << "Error: This in not a model correcting the star position (" << c << ")" << endl;
158 return;
159 }
160
161 Clear();
162
163 cout << endl;
164
165 Double_t val;
166 fin >> val;
167 cout << "Number of observed stars: " << val << endl;
168 fin >> val;
169 cout << "Sky RMS: " << val << "\"" << endl;
170 fin >> val;
171 cout << "Refraction Constant A: " << val << "\"" << endl;
172 fin >> val;
173 cout << "Refraction Constant B: " << val << "\"" << endl;
174
175 cout << endl;
176
177 cout << " & = Name Value Sigma" << endl;
178 cout << "--------------------------------------------------" << endl;
179
180 while (fin)
181 {
182 TString str;
183 fin >> str;
184
185 if (str=="END")
186 break;
187
188 if (str[0]=='&')
189 {
190 cout << " & ";
191 str.Remove(0);
192 }
193 else
194 cout << " ";
195
196 if (str[1]=='=')
197 {
198 cout << "= ";
199 str.Remove(0);
200 }
201 else
202 cout << " ";
203
204 fin >> val;
205 cout << str << "\t" << setw(11) << val << "° \t";
206 val *= kDeg2Rad;
207
208 // Find parameter
209 Int_t n = -1;
210 for (int i=0; i<fNumPar; i++)
211 if (str==fName[i])
212 {
213 n = i;
214 *fCoeff[i] = val;
215 break;
216 }
217
218 fin >> val;
219 cout << setw(9) << val << "°" << endl;
220
221 // corresponding error
222 fError[n] = val*kDeg2Rad;
223 }
224 cout << endl;
225}
226
227void MBending::Save(const char *name)
228{
229 /*
230 ! MMT 1987 July 8
231 ! T 36 7.3622 41.448 -0.0481
232 ! IA -37.5465 20.80602
233 ! IE -13.9180 1.25217
234 ! NPAE +7.0751 26.44763
235 ! CA -6.9149 32.05358
236 ! AN +0.5053 1.40956
237 ! AW -2.2016 1.37480
238 ! END
239 */
240
241 ofstream fout(name);
242 if (!fout)
243 {
244 cout << "Error: Cannot open file '" << name << "'" << endl;
245 return;
246 }
247
248 MTime t;
249 t.Now();
250
251 fout << "MAGIC1 " << t << endl;
252 fout << "S 00 000000 000000 0000000" << endl;
253 fout << setprecision(8);
254 for (int i=0; i<fNumPar; i++)
255 {
256 fout << " " << setw(6) << GetName(i) << " ";
257 fout << setw(13) << *fCoeff[i]*kRad2Deg << " ";
258 fout << setw(11) << fError[i]*kRad2Deg << endl;
259 }
260 fout << "END" << endl;
261}
262
263Double_t MBending::Sign(Double_t val, Double_t alt)
264{
265 // Some pointing corrections are defined as Delta ZA, which
266 // is (P. Wallace) defined [0,90]deg while Alt is defined
267 // [0,180]deg
268 return (TMath::Pi()/2-alt < 0 ? -val : val);
269}
270
271AltAz MBending::AddOffsets(const AltAz &aa) const
272{
273 // Correct [rad]
274 // zdaz [rad]
275 AltAz p = aa;
276
277 const AltAz I(fIe, fIa);
278 p += I;
279
280 return p;
281}
282
283AltAz MBending::SubtractOffsets(const AltAz &aa) const
284{
285 // Correct [rad]
286 // zdaz [rad]
287 AltAz p = aa;
288
289 const AltAz I(fIe, fIa);
290 p -= I;
291
292 return p;
293}
294
295AltAz MBending::CalcAnAw(const AltAz &p, Int_t sign) const
296{
297 // Corrections for AN and AW without approximations
298 // as done by Patrick Wallace. The approximation cannot
299 // be used for MAGIC because the correctioon angle
300 // AW (~1.5deg) is not small enough.
301
302 // Vector in cartesian coordinates
303 TVector3 v1;
304
305 // Set Azimuth and Elevation
306 v1.SetMagThetaPhi(1, TMath::Pi()/2-p.Alt(), p.Az());
307
308
309 TVector3 v2(v1);
310// cout << sign << endl;
311
312// cout << "v1: " << v1.Theta()*TMath::RadToDeg() << " " << v1.Phi()*TMath::RadToDeg() << endl;
313
314 // Rotate around the x- and y-axis
315 v1.RotateY(sign*fAn);
316 v1.RotateX(sign*fAw);
317
318// cout << "v1: " << v1.Theta()*TMath::RadToDeg() << " " << v1.Phi()*TMath::RadToDeg() << endl;
319// cout << "v2: " << v2.Theta()*TMath::RadToDeg() << " " << v2.Theta()*TMath::RadToDeg() << endl;
320
321 // cout << "dv: " << (v2.Theta()-v1.Theta())*TMath::RadToDeg() << " " << (v2.Phi()-v1.Phi())*TMath::RadToDeg() << endl;
322
323 Double_t dalt = v1.Theta()-v2.Theta();
324 Double_t daz = v1.Phi() -v2.Phi();
325
326 //cout << dalt*TMath::RadToDeg() << " " << daz*TMath::RadToDeg() << endl;
327
328 if (daz>TMath::Pi())
329 daz -= TMath::TwoPi();
330 if (daz<-TMath::Pi())
331 daz += TMath::TwoPi();
332
333// if (daz>TMath::Pi()/2)
334// {
335// }
336
337 AltAz d(dalt, daz);
338 return d;
339
340 // Calculate Delta Azimuth and Delta Elevation
341 /*
342 AltAz d(TMath::Pi()/2-v1.Theta(), v1.Phi());
343
344 cout << "p : " << p.Alt()*TMath::RadToDeg() << " " << p.Az()*TMath::RadToDeg() << endl;
345 cout << "d : " << d.Alt()*TMath::RadToDeg() << " " << d.Az()*TMath::RadToDeg() << endl;
346 d -= p;
347 cout << "d-p: " << d.Alt()*TMath::RadToDeg() << " " << d.Az()*TMath::RadToDeg() << endl;
348 d *= sign;
349 cout << "d* : " << d.Alt()*TMath::RadToDeg() << " " << d.Az()*TMath::RadToDeg() << endl;
350
351
352 cout << "p2: " << 90-p.Alt()*TMath::RadToDeg() << " " << p.Az()*TMath::RadToDeg() << endl;
353 cout << "d2: " << 90-d.Alt()*TMath::RadToDeg() << " " << d.Az()*TMath::RadToDeg() << endl;
354
355 Int_t s1 = 90-d.Alt()*TMath::RadToDeg() < 0 ? -1 : 1;
356 Int_t s2 = 90-p.Alt()*TMath::RadToDeg() < 0 ? -1 : 1;
357
358
359 if (s1 != s2)
360 {
361 //90-d.Alt() <-- -90+d.Alt()
362
363 d.Alt(d.Alt()-TMath::Pi());
364 cout << "Alt-" << endl;
365 }
366 cout << "d': " << 90-d.Alt()*TMath::RadToDeg() << " " << d.Az()*TMath::RadToDeg() << endl;
367 /*
368 // Fix 'direction' of output depending on input vector
369 if (TMath::Pi()/2-sign*p.Alt()<0)
370 {
371 d.Alt(d.Alt()-TMath::Pi());
372 cout << "Alt-" << endl;
373 }
374 //if (TMath::Pi()/2-sign*p.Alt()>TMath::Pi())
375 //{
376 // d.Alt(TMath::Pi()-d.Alt());
377 // cout << "Alt+" << endl;
378 //}
379
380 // Align correction into [-180,180]
381 while (d.Az()>TMath::Pi())
382 {
383 d.Az(d.Az()-TMath::Pi()*2);
384 cout << "Az-" << endl;
385 }
386 while (d.Az()<-TMath::Pi())
387 {
388 d.Az(d.Az()+TMath::Pi()*2);
389 cout << "Az+" << endl;
390 }
391 */
392 return d;
393}
394
395
396AltAz MBending::Correct(const AltAz &aa) const
397{
398 // Correct [rad]
399 // zdaz [rad]
400 AltAz p = aa;
401
402 DEBUG(cout << setprecision(16));
403 DEBUG(cout << "Bend7: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
404
405 const AltAz CRX(-fCrx*sin(p.Az()-p.Alt()), fCrx*cos(p.Az()-p.Alt())/cos(p.Alt()));
406 const AltAz CRY(-fCry*cos(p.Az()-p.Alt()), -fCry*sin(p.Az()-p.Alt())/cos(p.Alt()));
407 p += CRX;
408 p += CRY;
409
410 DEBUG(cout << "Bend6: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
411
412 const AltAz NRX(fNrx*sin(p.Alt()), -fNrx);
413 const AltAz NRY(fNry*cos(p.Alt()), -fNry*tan(p.Alt()));
414 p += NRX;
415 p += NRY;
416
417 DEBUG(cout << "Bend5: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
418
419 const AltAz CES(-fEces*sin(p.Alt()), -fAces*sin(p.Az()));
420 const AltAz CEC(-fEcec*cos(p.Alt()), -fAcec*cos(p.Az()));
421 p += CES;
422 p += CEC;
423
424 DEBUG(cout << "Bend4: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
425
426 const AltAz TX(Sign(fTx/tan(p.Alt()), p.Alt()), 0);
427 const AltAz TF(Sign(fTf*cos(p.Alt()), p.Alt()), 0);
428 //p += TX;
429 p += TF;
430
431
432 DEBUG(cout << "Bend3: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
433
434 /*
435 //New Corrections for NPAE and CA:
436 TVector3 v(1.,1.,1.); // Vector in cartesian coordinates
437
438 //Set Azimuth and Elevation
439 v.SetPhi(p.Az());
440 v.SetTheta(TMath::Pi()/2-p.Alt());
441 //Rotation Vectors:
442 TVector3 vNpae( cos(p.Az()), sin(p.Az()), 0);
443 TVector3 vCa( -cos(p.Az())*cos(p.Alt()), -sin(p.Az())*cos(p.Alt()), sin(p.Alt()));
444 //Rotate around the vectors vNpae and vCa
445 v.Rotate(fNpae, vNpae);
446 v.Rotate(fCa, vCa);
447
448 p.Az(v.Phi());
449 p.Alt(TMath::Pi()/2-v.Theta());
450 */
451
452 //Old correction terms for Npae and Ca:
453 const AltAz CA(0, -fCa/cos(p.Alt()));
454 p += CA;
455
456 const AltAz NPAE(0, -fNpae*tan(p.Alt()));
457 p += NPAE;
458
459 DEBUG(cout << "Bend2: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
460
461 const AltAz ANAW(CalcAnAw(p, -1));
462 p += ANAW;
463
464 /* Old correction terms for An and Aw:
465 const AltAz AW( fAw*sin(p.Az()), -fAw*cos(p.Az())*tan(p.Alt()));
466 const AltAz AN(-fAn*cos(p.Az()), -fAn*sin(p.Az())*tan(p.Alt()));
467 p += AW;
468 p += AN;
469 */
470
471 DEBUG(cout << "Bend1: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
472
473 const AltAz FLOP(Sign(fFlop, p.Alt()), 0);
474 p += FLOP;
475
476 const AltAz MAGIC1(fMagic1*TMath::Sign(1., sin(p.Az())), 0);
477 p += MAGIC1;
478
479 const AltAz I(fIe, fIa);
480 p += I;
481
482 DEBUG(cout << "Bend0: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
483
484 return p;
485}
486
487AltAz MBending::CorrectBack(const AltAz &aa) const
488{
489 // Correct [rad]
490 // zdaz [rad]
491 AltAz p = aa;
492
493 DEBUG(cout << setprecision(16));
494 DEBUG(cout << "Back0: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
495
496 const AltAz I(fIe, fIa);
497 p -= I;
498
499 const AltAz MAGIC1(fMagic1*TMath::Sign(1., sin(p.Az())), 0);
500 p -= MAGIC1;
501
502 //const AltAz MAGIC1(fMagic1*sin(p.Az()), 0);
503 //p -= MAGIC1;
504
505 const AltAz FLOP(Sign(fFlop, p.Alt()), 0);
506 p -= FLOP;
507
508 DEBUG(cout << "Back1: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
509
510 /* Old correction terms for An and Aw:
511 const AltAz AN(-fAn*cos(p.Az()), -fAn*sin(p.Az())*tan(p.Alt()));
512 const AltAz AW( fAw*sin(p.Az()), -fAw*cos(p.Az())*tan(p.Alt()));
513 p -= AN;
514 p -= AW;
515 */
516
517 const AltAz ANAW(CalcAnAw(p, -1));
518 p -= ANAW;
519
520 DEBUG(cout << "Back2: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
521
522 //Old Correction terms for Npae and Ca:
523 const AltAz NPAE(0, -fNpae*tan(p.Alt()));
524 p -= NPAE;
525
526 const AltAz CA(0, -fCa/cos(p.Alt()));
527 p -= CA;
528
529 /*
530 //New Correction term for Npae and Ca:
531 TVector3 v2(1.,1.,1.); // Vector in cartesian coordinates
532 //Set Azimuth and Elevation
533 v2.SetPhi(p.Az());
534 v2.SetTheta(TMath::Pi()/2-p.Alt());
535 //Rotation Vectors:
536 TVector3 vNpae( cos(p.Az()), sin(p.Az()), 0);
537 TVector3 vCa( -cos(p.Az())*cos(p.Alt()), -sin(p.Az())*cos(p.Alt()), sin(p.Alt()));
538 //Rotate around the vectors vCa and vNpae
539 v2.Rotate(-fCa, vCa);
540 v2.Rotate(-fNpae, vNpae);
541
542 p.Az(v2.Phi());
543 p.Alt(TMath::Pi()/2-v2.Theta());
544 */
545
546 DEBUG(cout << "Back3: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
547
548 const AltAz TF(Sign(fTf*cos(p.Alt()), p.Alt()), 0);
549 const AltAz TX(Sign(fTx/tan(p.Alt()), p.Alt()), 0);
550 p -= TF;
551 //p -= TX;
552
553 DEBUG(cout << "Back4: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
554
555 const AltAz CEC(-fEcec*cos(p.Alt()), -fAcec*cos(p.Az()));
556 const AltAz CES(-fEces*sin(p.Alt()), -fAces*sin(p.Az()));
557 p -= CEC;
558 p -= CES;
559
560 DEBUG(cout << "Back5: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
561
562 const AltAz NRY(fNry*cos(p.Alt()), -fNry*tan(p.Alt()));
563 const AltAz NRX(fNrx*sin(p.Alt()), -fNrx);
564 p -= NRY;
565 p -= NRX;
566
567 DEBUG(cout << "Back6: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
568
569 const AltAz CRY(-fCry*cos(p.Az()-p.Alt()), -fCry*sin(p.Az()-p.Alt())/cos(p.Alt()));
570 const AltAz CRX(-fCrx*sin(p.Az()-p.Alt()), fCrx*cos(p.Az()-p.Alt())/cos(p.Alt()));
571 p -= CRY;
572 p -= CRX;
573
574 DEBUG(cout << "Back7: " << 90-p.Alt()*180/TMath::Pi() << " " << p.Az()*180/TMath::Pi() << endl);
575
576 return p;
577}
578
579ZdAz MBending::Correct(const ZdAz &zdaz) const
580{
581 // Correct [rad]
582 // zdaz [rad]
583 AltAz p(TMath::Pi()/2-zdaz.Zd(), zdaz.Az());
584 AltAz c = Correct(p);
585 return ZdAz(TMath::Pi()/2-c.Alt(), c.Az());
586}
587
588ZdAz MBending::CorrectBack(const ZdAz &zdaz) const
589{
590 // Correct [rad]
591 // zdaz [rad]
592 AltAz p(TMath::Pi()/2-zdaz.Zd(), zdaz.Az());
593 AltAz c = CorrectBack(p);
594 return ZdAz(TMath::Pi()/2-c.Alt(), c.Az());
595}
596
597void MBending::SetParameters(const Double_t *par, Int_t n)
598{
599 Clear();
600
601 while (n--)
602 *fCoeff[n] = par[n]/kRad2Deg;
603}
604
605void MBending::GetParameters(Double_t *par, Int_t n) const
606{
607 while (n--)
608 par[n] = *fCoeff[n]*kRad2Deg;
609}
610
611void MBending::SetMinuitParameters(TMinuit &m, Int_t n) const
612{
613 if (n<0)
614 n = fNumPar;
615
616 Int_t ierflg = 0;
617
618 while (n--)
619 m.mnparm(n, fName[n], *fCoeff[n]*kRad2Deg, 1, -360, 360, ierflg);
620}
621
622void MBending::GetMinuitParameters(TMinuit &m, Int_t n)
623{
624 if (n<0 || n>m.GetNumPars())
625 n = m.GetNumPars();
626
627 while (n--)
628 {
629 m.GetParameter(n, *fCoeff[n], fError[n]);
630 *fCoeff[n] /= kRad2Deg;
631 fError[n] /= kRad2Deg;
632 }
633}
634/*
635void FormatPar(TMinuit &m, Int_t n)
636{
637 Double_t par, err;
638 m.GetParameter(n, par, err);
639
640 int expp = (int)log10(par);
641 int expe = (int)log10(err);
642
643 if (err<2*pow(10, expe))
644 expe--;
645
646 Int_t exp = expe>expp ? expp : expe;
647
648 par = (int)(par/pow(10, exp)) * pow(10, exp);
649 err = (int)(err/pow(10, exp)) * pow(10, exp);
650
651 cout << par << " +- " << err << flush;
652}
653*/
654void MBending::PrintMinuitParameters(TMinuit &m, Int_t n) const
655{
656 if (n<0)
657 n = m.GetNumPars();
658
659 cout << setprecision(3);
660
661 Double_t par, err;
662
663 while (n--)
664 {
665 m.GetParameter(n, par, err);
666 cout << Form(" %2d %6s: ", n, (const char*)fName[n]);
667 cout << setw(8) << par << " \xb1 " << setw(6) << err << endl;
668 }
669}
Note: See TracBrowser for help on using the repository browser.