source: trunk/MagicSoft/Cosy/main/MBending.cc@ 6724

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