source: trunk/MagicSoft/Cosy/tpoint/gui.C@ 2977

Last change on this file since 2977 was 2568, checked in by tbretz, 21 years ago
*** empty log message ***
File size: 25.9 KB
Line 
1#include <fstream.h>
2#include <iostream.h>
3#include <iomanip.h>
4
5#include <TGFrame.h>
6#include <TGLabel.h>
7#include <TGButton.h>
8
9#include <TF1.h>
10#include <TH1.h>
11#include <TH2.h>
12#include <TGraph.h>
13
14#include <TList.h>
15#include <TStyle.h>
16#include <TMinuit.h>
17
18#include <TView.h>
19#include <TLine.h>
20#include <TMarker.h>
21#include <TCanvas.h>
22
23#include "coord.h"
24
25#include "MGList.h"
26#include "MBending.h"
27
28class Set : public TObject
29{
30 friend ifstream &operator>>(ifstream &fin, Set &set);
31private:
32 Double_t fStarAz;
33 Double_t fStarEl;
34
35 Double_t fRawAz;
36 Double_t fRawEl;
37
38public:
39 Set(Double_t sel=0, Double_t saz=0, Double_t rel=0, Double_t raz=0) :
40 fStarAz(saz*TMath::Pi()/180),
41 fStarEl(sel*TMath::Pi()/180),
42 fRawAz(raz*TMath::Pi()/180),
43 fRawEl(rel*TMath::Pi()/180)
44 {
45 }
46
47 Double_t GetResidual() const
48 {
49 Double_t del = fRawEl-fStarEl;
50 Double_t daz = fRawAz-fStarAz;
51 Double_t dphi2 = daz/2.;
52 Double_t cos2 = cos(dphi2)*cos(dphi2);
53 Double_t sin2 = sin(dphi2)*sin(dphi2);
54 Double_t d = cos(del)*cos2 - cos(fRawEl+fStarEl)*sin2;
55
56 Double_t dist = acos(d);
57
58 return dist * 180 / TMath::Pi();
59 }
60
61 void operator=(Set &set)
62 {
63 fStarAz = set.fStarAz;
64 fStarEl = set.fStarEl;
65 fRawAz = set.fRawAz;
66 fRawEl = set.fRawEl;
67 }
68
69 Double_t GetDEl() const { return (fRawEl-fStarEl)*180/TMath::Pi(); }
70 Double_t GetDZd() const { return -GetDEl(); }
71 Double_t GetDAz() const { return (fRawAz-fStarAz)*180/TMath::Pi(); }
72 Double_t GetStarEl() const { return fStarEl*180/TMath::Pi(); }
73 Double_t GetStarZd() const { return 90.-fStarEl*180/TMath::Pi(); }
74 Double_t GetStarAz() const { return fStarAz*180/TMath::Pi(); }
75 Double_t GetRawEl() const { return fRawEl*180/TMath::Pi(); }
76 Double_t GetRawAz() const { return fRawAz*180/TMath::Pi(); }
77 Double_t GetRawZd() const { return 90.-fRawEl*180/TMath::Pi(); }
78
79 ZdAz GetStarZdAz() const { return ZdAz(TMath::Pi()/2-fStarEl, fStarAz); }
80 AltAz GetStarAltAz() const { return AltAz(fStarEl, fStarAz); }
81
82 ZdAz GetRawZdAz() const { return ZdAz(TMath::Pi()/2-fRawEl, fRawAz); }
83 AltAz GetRawAltAz() const { return AltAz(fRawEl, fRawAz); }
84
85 void AdjustEl(Double_t del) { fStarEl += del*TMath::Pi()/180; }
86 void AdjustAz(Double_t daz) { fStarAz += daz*TMath::Pi()/180; }
87
88 void Adjust(const MBending &bend)
89 {
90 AltAz p = bend(GetStarAltAz());
91 fStarEl = p.Alt();
92 fStarAz = p.Az();
93 }
94 void AdjustBack(const MBending &bend)
95 {
96 AltAz p = bend.CorrectBack(GetRawAltAz());
97 fRawEl = p.Alt();
98 fRawAz = p.Az();
99 }
100 ClassDef(Set, 0)
101};
102
103ClassImp(Set);
104
105ifstream &operator>>(ifstream &fin, Set &set)
106{
107 Double_t v[4];
108 fin >> v[0];
109 fin >> v[1];
110 fin >> v[2];
111 fin >> v[3];
112
113 Double_t dummy;
114 fin>>dummy;
115 fin>>dummy;
116 fin>>dummy;
117
118 set.fStarAz = v[0]*TMath::Pi()/180;
119 set.fStarEl = v[1]*TMath::Pi()/180;
120
121 set.fRawAz = v[2]*TMath::Pi()/180;
122 set.fRawEl = v[3]*TMath::Pi()/180;
123
124 return fin;
125}
126
127class MFit : public TGMainFrame
128{
129private:
130 enum {
131 kTbFit = 19, //MBending::GetNumPar(), // FIXME!!!
132 kTbLoad,
133 kTbSave,
134 kTbLoadStars,
135 kTbReset,
136 kTbResetStars
137 };
138
139 MGList *fList;
140
141 TList fOriginal;
142 TList fCoordinates;
143 TList fLabel;
144
145 MBending fBending;
146
147 FontStruct_t fFont;
148
149 void Fcn(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
150 {
151 f = 0;
152
153 MBending bend;
154 bend.SetParameters(par); // Set Parameters [deg] to MBending
155
156 for (int i=0; i<fCoordinates.GetSize(); i++)
157 {
158 Set set = *(Set*)fCoordinates.At(i);
159
160 set.Adjust(bend);
161
162 Double_t err = 0.02; // [deg] = 1SE
163 Double_t res = set.GetResidual()/err;
164
165 f += res*res;
166 }
167
168 //f /= (fCoordinates.GetSize()-7)*(fCoordinates.GetSize()-7);
169 //f /= fCoordinates.GetSize()*fCoordinates.GetSize();
170 f /= fCoordinates.GetSize();
171 }
172
173 static void fcn(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
174 {
175 ((MFit*)gMinuit->GetObjectFit())->Fcn(npar, gin, f, par, iflag);
176 }
177
178 void DrawMarker(TVirtualPad *pad, Double_t r0, Double_t phi0, Double_t r1, Double_t phi1)
179 {
180 TView *view = pad->GetView();
181
182 if (!view)
183 {
184 cout << "No View!" << endl;
185 return;
186 }
187
188 TMarker mark0;
189 TMarker mark1;
190 mark0.SetMarkerStyle(kStar);
191 mark1.SetMarkerStyle(kStar);
192 mark1.SetMarkerColor(kRed);
193
194 r0 /= 90;
195 r1 /= 90;
196 phi0 *= TMath::Pi()/180;
197 phi1 *= TMath::Pi()/180;
198
199 Double_t x0[3] = { r0*cos(phi0), r0*sin(phi0), 0};
200 Double_t x1[3] = { r1*cos(phi1), r1*sin(phi1), 0};
201
202 Double_t y0[3], y1[3];
203
204 view->WCtoNDC(x0, y0);
205 view->WCtoNDC(x1, y1);
206
207 mark0.DrawMarker(y0[0], y0[1]);
208 mark1.DrawMarker(y1[0], y1[1]);
209 }
210
211 void DrawPolLine(TVirtualPad *pad, Double_t r0, Double_t phi0, Double_t r1, Double_t phi1)
212 {
213 TView *view = pad->GetView();
214
215 if (!view)
216 {
217 cout << "No View!" << endl;
218 return;
219 }
220 /*
221 if (r0<0)
222 {
223 r0 = -r0;
224 phi0 += 180;
225 }
226 if (r1<0)
227 {
228 r1 = -r1;
229 phi1 += 180;
230 }
231
232 phi0 = fmod(phi0+360, 360);
233 phi1 = fmod(phi1+360, 360);
234
235 if (phi1-phi0<-180)
236 phi1+=360;
237 */
238 TLine line;
239 line.SetLineWidth(2);
240 line.SetLineColor(kBlue);
241
242 Double_t p0 = phi0<phi1?phi0:phi1;
243 Double_t p1 = phi0<phi1?phi1:phi0;
244
245 if (phi0>phi1)
246 {
247 Double_t d = r1;
248 r1 = r0;
249 r0 = d;
250 }
251
252 r0 /= 90;
253 r1 /= 90;
254
255 Double_t dr = r1-r0;
256 Double_t dp = p1-p0;
257
258 Double_t x0[3] = { r0*cos(p0*TMath::Pi()/180), r0*sin(p0*TMath::Pi()/180), 0};
259
260 for (double i=p0+10; i<p1+10; i+=10)
261 {
262 if (i>p1)
263 i=p1;
264
265 Double_t r = dr/dp*(i-p0)+r0;
266 Double_t p = TMath::Pi()/180*i;
267
268 Double_t x1[3] = { r*cos(p), r*sin(p), 0};
269
270 Double_t y0[3], y1[3];
271
272 view->WCtoNDC(x0, y0);
273 view->WCtoNDC(x1, y1);
274
275 line.DrawLine(y0[0], y0[1], y1[0], y1[1]);
276
277 x0[0] = x1[0];
278 x0[1] = x1[1];
279 }
280 }
281
282 void DrawSet(TVirtualPad *pad, Set &set, Float_t scale=-1)
283 {
284 Double_t r0 = set.GetRawZd();
285 Double_t phi0 = set.GetRawAz();
286 Double_t r1 = set.GetStarZd();
287 Double_t phi1 = set.GetStarAz();
288
289 if (r0<0)
290 {
291 r0 = -r0;
292 phi0 += 180;
293 }
294 if (r1<0)
295 {
296 r1 = -r1;
297 phi1 += 180;
298 }
299
300 phi0 = fmod(phi0+360, 360);
301 phi1 = fmod(phi1+360, 360);
302
303 if (phi1-phi0<-180)
304 phi1+=360;
305
306 if (scale<0 || scale>1000)
307 scale = -1;
308
309 if (scale>0)
310 {
311 Double_t d = r1-r0;
312 r0 += scale*d;
313 r1 -= scale*d;
314 d = phi1-phi0;
315 phi0 += scale*d;
316 phi1 -= scale*d;
317
318 DrawPolLine(pad, r0, phi0, r1, phi1);
319 DrawMarker(pad, r0, phi0, r1, phi1);
320 }
321 else
322 DrawMarker(pad, r1, phi1, 0 ,0);
323 }
324
325 void Fit()
326 {
327 if (fOriginal.GetSize()==0)
328 {
329 cout << "Sorry, no input data loaded..." << endl;
330 return;
331 }
332
333 fCoordinates.Delete();
334 for (int i=0; i<fOriginal.GetSize(); i++)
335 fCoordinates.Add(new Set(*(Set*)fOriginal.At(i)));
336
337 cout << "-----------------------------------------------------------------------" << endl;
338
339 TH1F hres1("Res1", " Residuals before correction ", fOriginal.GetSize()/2, 0, 3);
340 TH1F hres2("Res2", " Residuals after correction ", fOriginal.GetSize(), 0, 3);
341 TH1F hres3("Res3", " Residuals after backward correction ", fOriginal.GetSize(), 0, 3);
342
343 hres1.SetXTitle("\\Delta [\\circ]");
344 hres1.SetYTitle("Counts");
345
346 hres2.SetXTitle("\\Delta [\\circ]");
347 hres2.SetYTitle("Counts");
348
349 hres3.SetXTitle("\\Delta [\\circ]");
350 hres3.SetYTitle("Counts");
351
352 TGraph gdaz;
353 TGraph gdzd;
354 TGraph gaz;
355 TGraph gzd;
356 TGraph graz;
357 TGraph grzd;
358
359 gdaz.SetTitle(" \\Delta Az vs. Zd ");
360 gdzd.SetTitle(" \\Delta Zd vs. Az ");
361
362 gaz.SetTitle(" \\Delta Az vs. Az ");
363 gzd.SetTitle(" \\Delta Zd vs. Zd ");
364
365 graz.SetTitle(" \\Delta vs. Az ");
366 grzd.SetTitle(" \\Delta vs. Zd ");
367
368 TMinuit minuit(MBending::GetNumPar()); //initialize TMinuit with a maximum of 5 params
369 minuit.SetObjectFit(this);
370 minuit.SetPrintLevel(-1);
371 minuit.SetFCN(fcn);
372
373 fBending.SetMinuitParameters(minuit, MBending::GetNumPar()); // Init Parameters [deg]
374
375 for (int i=0; i<MBending::GetNumPar(); i++)
376 {
377 TGButton *l = (TGButton*)fList->FindWidget(i);
378 minuit.FixParameter(i);
379 if (l->GetState()==kButtonDown)
380 minuit.Release(i);
381 }
382
383 //minuit.Command("SHOW PARAMETERS");
384 //minuit.Command("SHOW LIMITS");
385
386 cout << endl;
387 cout << "Starting fit..." << endl;
388 cout << "For the fit an measurement error in the residual of ";
389 cout << "0.02deg (=1SE) is assumed." << endl;
390 cout << endl;
391
392 Int_t ierflg = 0;
393 ierflg = minuit.Migrad();
394 cout << "Migrad returns " << ierflg << endl;
395 // minuit.Release(2);
396 ierflg = minuit.Migrad();
397 cout << "Migrad returns " << ierflg << endl << endl;
398
399 //
400 // Get Fit Results
401 //
402 fBending.GetMinuitParameters(minuit);
403 fBending.PrintMinuitParameters(minuit);
404 //fBending.Save("bending_magic.txt");
405
406
407 //
408 // Make a copy of all list entries
409 //
410 TList list;
411 list.SetOwner();
412 for (int i=0; i<fCoordinates.GetSize(); i++)
413 list.Add(new Set(*(Set*)fCoordinates.At(i)));
414
415 //
416 // Correct for Offsets only
417 //
418 TArrayD par;
419 fBending.GetParameters(par);
420 for (int i=2; i<MBending::GetNumPar(); i++)
421 par[i]=0;
422
423 MBending b2;
424 b2.SetParameters(par);
425
426 //
427 // Calculate correction and residuals
428 //
429 for (int i=0; i<fCoordinates.GetSize(); i++)
430 {
431 Set &set0 = *(Set*)fCoordinates.At(i);
432
433 ZdAz za = set0.GetStarZdAz()*kRad2Deg;
434
435 //
436 // Correct for offsets only
437 //
438 Set set1(set0);
439 set1.Adjust(b2);
440
441 hres1.Fill(set1.GetResidual());
442
443 set0.Adjust(fBending);
444 hres2.Fill(set0.GetResidual());
445
446 Double_t dz = fmod(set0.GetDAz()+720, 360);
447 if (dz>180)
448 dz -= 360;
449
450 gdzd.SetPoint(i, za.Az(), set0.GetDZd());
451 gdaz.SetPoint(i, za.Zd(), dz);
452 graz.SetPoint(i, za.Az(), set0.GetResidual());
453 grzd.SetPoint(i, za.Zd(), set0.GetResidual());
454 gaz.SetPoint( i, za.Az(), dz);
455 gzd.SetPoint( i, za.Zd(), set0.GetDZd());
456
457 }
458
459 //
460 // Check for overflows
461 //
462 const Stat_t ov = hres2.GetBinContent(hres2.GetNbinsX()+1);
463 if (ov>0)
464 cout << "WARNING: " << ov << " overflows in residuals." << endl;
465
466
467
468 cout << dec << endl;
469 cout << " Number of calls to FCN: " << minuit.fNfcn << endl;
470 cout << "Minimum value found for FCN (Chi^2?): " << minuit.fAmin << endl;
471 cout << " Fit-Probability(?): " << TMath::Prob(minuit.fAmin/*fOriginal.GetSize()*/, fOriginal.GetSize()-minuit.GetNumFreePars())*100 << "%" << endl;
472 cout << " Chi^2/NDF: " << minuit.fAmin/(fOriginal.GetSize()-minuit.GetNumFreePars()) << endl;
473 //cout << "Prob(?): " << TMath::Prob(fChisquare,ndf);
474
475
476
477 //
478 // Print all data sets for which the backward correction is
479 // twice times worse than the residual gotten from the
480 // bending correction itself
481 //
482 cout << endl;
483 cout << "Checking backward correction (raw-->star):" << endl;
484 for (int i=0; i<fCoordinates.GetSize(); i++)
485 {
486 Set set0(*(Set*)list.At(i));
487 Set &set1 = *(Set*)list.At(i);
488
489 set0.AdjustBack(fBending);
490 set1.Adjust(fBending);
491
492 const Double_t res0 = set0.GetResidual();
493 const Double_t res1 = set1.GetResidual();
494 const Double_t diff = TMath::Abs(res0-res1);
495
496 hres3.Fill(res0);
497
498 if (diff<hres2.GetMean()*0.66)
499 continue;
500
501 cout << "DBack: " << setw(6) << set0.GetStarZd() << " " << setw(7) << set0.GetStarAz() << ": ";
502 cout << "ResB="<< setw(7) << res0*60 << " ResF=" << setw(7) << res1*60 << " |ResB-ResF|=" << setw(7) << diff*60 << " arcmin" << endl;
503 }
504 cout << "OK." << endl;
505 cout << endl;
506
507
508
509
510 TCanvas *c1;
511
512 if ((c1 = (TCanvas*)gROOT->FindObject("CanvGraphs")))
513 delete c1;
514 if ((c1 = (TCanvas*)gROOT->FindObject("CanvResiduals")))
515 delete c1;
516
517
518
519 c1=new TCanvas("CanvGraphs", "Graphs");
520 c1->Divide(2,3,0,0);
521
522 c1->cd(1);
523 gPad->SetBorderMode(0);
524 TGraph *g=(TGraph*)gaz.DrawClone("A*");
525 g->SetBit(kCanDelete);
526 g->GetHistogram()->SetXTitle("Az [\\circ]");
527 g->GetHistogram()->SetYTitle("\\Delta Az [\\circ]");
528
529 c1->cd(2);
530 gPad->SetBorderMode(0);
531 g=(TGraph*)gdaz.DrawClone("A*");
532 g->SetBit(kCanDelete);
533 g->GetHistogram()->SetXTitle("Zd [\\circ]");
534 g->GetHistogram()->SetYTitle("\\Delta Az [\\circ]");
535 cout << "Mean dAz: " << g->GetMean(2) << " \xb1 " << g->GetRMS(2) << endl;
536
537 c1->cd(3);
538 gPad->SetBorderMode(0);
539 g=(TGraph*)gdzd.DrawClone("A*");
540 g->SetBit(kCanDelete);
541 g->GetHistogram()->SetXTitle("Az [\\circ]");
542 g->GetHistogram()->SetYTitle("\\Delta Zd [\\circ]");
543 cout << "Mean dZd: " << g->GetMean(2) << " \xb1 " << g->GetRMS(2) << endl;
544 cout << endl;
545
546 c1->cd(4);
547 gPad->SetBorderMode(0);
548 g=(TGraph*)gzd.DrawClone("A*");
549 g->SetBit(kCanDelete);
550 g->GetHistogram()->SetXTitle("Zd [\\circ]");
551 g->GetHistogram()->SetYTitle("\\Delta Zd [\\circ]");
552
553 c1->cd(5);
554 gPad->SetBorderMode(0);
555 g=(TGraph*)graz.DrawClone("A*");
556 g->SetBit(kCanDelete);
557 g->GetHistogram()->SetXTitle("Az [\\circ]");
558 g->GetHistogram()->SetYTitle("\\Delta [\\circ]");
559
560 c1->cd(6);
561 gPad->SetBorderMode(0);
562 g=(TGraph*)grzd.DrawClone("A*");
563 g->SetBit(kCanDelete);
564 g->GetHistogram()->SetXTitle("Zd [\\circ]");
565 g->GetHistogram()->SetYTitle("\\Delta [\\circ]");
566
567
568
569 //
570 // Print out the residual before and after correction in several
571 // units
572 //
573 cout << fCoordinates.GetSize() << " data sets." << endl << endl;
574 cout << "Total Spread of Residual:" << endl;
575 cout << "-------------------------" << endl;
576 cout << "before: " << hres1.GetMean() << " \xb1 " << hres1.GetRMS() << " deg " << endl;
577 cout << "after: " << hres2.GetMean() << " \xb1 " << hres2.GetRMS() << " deg " << endl;
578 cout << "before: " << (int)(hres1.GetMean()*60+.5) << " \xb1 " << (int)(hres1.GetRMS()*60+.5) << " arcmin" << endl;
579 cout << "after: " << (int)(hres2.GetMean()*60+.5) << " \xb1 " << (int)(hres2.GetRMS()*60+.5) << " arcmin" << endl;
580 cout << "before: " << (int)(hres1.GetMean()*60*60/23.4+.5) << " \xb1 " << (int)(hres1.GetRMS()*60*60/23.4+.5) << " pix" << endl;
581 cout << "after: " << (int)(hres2.GetMean()*60*60/23.4+.5) << " \xb1 " << (int)(hres2.GetRMS()*60*60/23.4+.5) << " pix" << endl;
582 cout << "after: " << (int)(hres2.GetMean()*16384/360+.5) << " \xb1 " << (int)(hres2.GetRMS()*16384/360+.5) << " SE" << endl;
583 cout << endl;
584 cout << "backw: " << (int)(hres3.GetMean()*60+.5) << " \xb1 " << (int)(hres3.GetRMS()*60+.5) << " arcmin" << endl;
585 cout << endl; // ±
586
587
588
589 gStyle->SetOptStat(1110);
590 gStyle->SetStatFormat("6.2g");
591
592
593 c1=new TCanvas("CanvResiduals", "Residuals", 800, 800);
594 c1->Divide(2, 2, 0, 0);
595
596 c1->cd(2);
597 gPad->SetBorderMode(0);
598 hres1.SetLineColor(kRed);
599 hres1.DrawCopy();
600
601 c1->cd(4);
602 gPad->SetBorderMode(0);
603 hres2.SetLineColor(kBlue);
604 TH1 *h=hres2.DrawCopy();
605 TF1 f("mygaus", "(gaus)", 0, 1);
606 f.SetLineColor(kMagenta/*6*/);
607 f.SetLineWidth(1);
608 f.SetParameter(0, h->GetBinContent(1));
609 f.FixParameter(1, 0);
610 f.SetParameter(2, h->GetRMS());
611 h->Fit("mygaus", "QR");
612 hres3.SetLineColor(kGreen);
613 hres3.SetLineStyle(kDashed);
614 hres3.DrawCopy("same");
615 cout << "Gaus-Fit Sigma: " << f.GetParameter(2) << "\xb0" << endl;
616 cout << "Fit-Probability: " << f.GetProb()*100 << "%" << endl;
617 cout << " Chi^2/NDF: " << f.GetChisquare()/f.GetNDF() << endl;
618
619 c1->cd(1);
620 gPad->SetBorderMode(0);
621 gPad->SetTheta(90);
622 gPad->SetPhi(90);
623 TH2F h2res1("Res2D1", " Dataset positions on the sky ", 32, 0, 360, 10, 0, 90);
624 h2res1.SetBit(TH1::kNoStats);
625 h2res1.DrawCopy("surf1pol");
626 gPad->Modified();
627 gPad->Update();
628 for (int i=0; i<fCoordinates.GetSize(); i++)
629 DrawSet(gPad, *(Set*)list.At(i));//, 10./hres1.GetMean());
630
631 c1->cd(3);
632 gPad->SetBorderMode(0);
633 gPad->SetTheta(90);
634 gPad->SetPhi(90);
635 h2res1.SetTitle(" Arb. Residuals after correction (scaled) ");
636 h2res1.DrawCopy("surf1pol");
637 gPad->Modified();
638 gPad->Update();
639 for (int i=0; i<fCoordinates.GetSize(); i++)
640 DrawSet(gPad, *(Set*)fCoordinates.At(i), 10./hres2.GetMean());
641
642 RaiseWindow();
643 }
644
645 void LoadStars()
646 {
647 const Int_t size = fOriginal.GetSize();
648
649 ifstream fin("tpoint.txt");
650
651 while (fin && fin.get()!='\n');
652 while (fin && fin.get()!='\n');
653 while (fin && fin.get()!='\n');
654 if (!fin)
655 {
656 cout << "File 'tpoint.txt' not found!" << endl;
657 return;
658 }
659
660 while (1)
661 {
662 Set set;
663
664 fin >> set; // Read data from file [deg], it is stored in [rad]
665
666 if (!fin)
667 break;
668
669 fOriginal.Add(new Set(set));
670 }
671
672 cout << "Found " << fOriginal.GetSize()-size << " sets of coordinates ";
673 cout << "(Total=" << fOriginal.GetSize() << ")" << endl;
674 }
675
676 Bool_t ProcessMessage(Long_t msg, Long_t mp1, Long_t mp2)
677 {
678 // cout << "Msg: " << hex << GET_MSG(msg) << endl;
679 // cout << "SubMsg: " << hex << GET_SUBMSG(msg) << dec << endl;
680 switch (GET_MSG(msg))
681 {
682 case kC_COMMAND:
683 switch (GET_SUBMSG(msg))
684 {
685 case kCM_BUTTON:
686 switch (mp1)
687 {
688 case kTbFit:
689 Fit();
690 DisplayBending();
691 return kTRUE;
692 case kTbLoad:
693 fBending.Load("bending_magic.txt");
694 DisplayBending();
695 return kTRUE;
696 case kTbSave:
697 fBending.Save("bending_magic.txt");
698 return kTRUE;
699 case kTbLoadStars:
700 LoadStars();
701 DisplayData();
702 return kTRUE;
703 case kTbReset:
704 fBending.Clear();
705 DisplayBending();
706 return kTRUE;
707 case kTbResetStars:
708 fOriginal.Delete();
709 DisplayData();
710 return kTRUE;
711 }
712 return kTRUE;
713 }
714 return kTRUE;
715 }
716 return kTRUE;
717 }
718
719 void AddTextButton(TGCompositeFrame *f, TString txt, Int_t id=-1, TGLayoutHints *h=0)
720 {
721 TGButton *but = new TGTextButton(f, txt, id);
722 but->Associate(this);
723 f->AddFrame(but, h);
724 fList->Add(but);
725
726 }
727
728 void AddCheckButton(TGCompositeFrame *f, TString txt, Int_t id=-1, TGLayoutHints *h=0)
729 {
730 TGButton *but = new TGCheckButton(f, txt, id);
731 but->Associate(this);
732 f->AddFrame(but, h);
733 fList->Add(but);
734 }
735
736 TGLabel *AddLabel(TGCompositeFrame *f, TString txt, TGLayoutHints *h=0)
737 {
738 TGLabel *l = new TGLabel(f, txt/*, TGLabel::GetDefaultGC()(), fFont*/);
739 f->AddFrame(l, h);
740 fList->Add(l);
741 fLabel.Add(l);
742 return l;
743 }
744
745 void DisplayBending()
746 {
747 TArrayD par, err;
748 fBending.GetParameters(par);
749 fBending.GetError(err);
750
751 TGLabel *l;
752
753 for (int i=0; i<MBending::GetNumPar(); i++)
754 {
755 l = (TGLabel*)fLabel.At(i);
756 l->SetText(Form("%.4f\xb0", par[i]));
757
758 l = (TGLabel*)fLabel.At(MBending::GetNumPar()+i);
759 l->SetText(Form("\xb1 %8.4f\xb0", err[i]));
760 }
761 }
762
763 void DisplayData()
764 {
765 TGLabel *l = (TGLabel*)fLabel.At(3*MBending::GetNumPar());
766 l->SetText(Form("%d data sets loaded.", fOriginal.GetSize()));
767 }
768
769public:
770 ~MFit()
771 {
772 if (fFont)
773 gVirtualX->DeleteFont(fFont);
774 }
775 MFit() : TGMainFrame(gClient->GetRoot(), 740, 370, kHorizontalFrame)
776 {
777 fCoordinates.SetOwner();
778 fOriginal.SetOwner();
779
780 fList = new MGList;
781 fList->SetOwner();
782
783 fFont = gVirtualX->LoadQueryFont("7x13bold");
784
785 TGLayoutHints *hints0 = new TGLayoutHints(kLHintsExpandY, 7, 5, 5, 6);
786 TGLayoutHints *hints1 = new TGLayoutHints(kLHintsExpandX|kLHintsExpandY, 5, 7, 5, 6);
787 TGLayoutHints *hints2 = new TGLayoutHints(kLHintsCenterX|kLHintsCenterY, 5, 5, 5, 5);
788 fList->Add(hints0);
789 fList->Add(hints1);
790 fList->Add(hints2);
791
792 TGGroupFrame *grp1 = new TGGroupFrame(this, "Control", kVerticalFrame);
793 AddFrame(grp1, hints0);
794 fList->Add(grp1);
795
796 TGGroupFrame *grp2 = new TGGroupFrame(this, "Parameters", kHorizontalFrame);
797 AddFrame(grp2, hints1);
798 fList->Add(grp2);
799
800
801
802 TGLayoutHints *hints4 = new TGLayoutHints(kLHintsExpandX, 5, 5, 5);
803 TGLayoutHints *hints5 = new TGLayoutHints(kLHintsExpandX, 5, 5, 15);
804 AddTextButton(grp1, "Load Pointing Model", kTbLoad, hints5);
805 AddTextButton(grp1, "Save Pointing Model", kTbSave, hints4);
806 AddTextButton(grp1, "Fit Parameters", kTbFit, hints5);
807 AddTextButton(grp1, "Reset Parameters", kTbReset, hints4);
808 AddTextButton(grp1, "Load Stars", kTbLoadStars, hints5);
809 AddTextButton(grp1, "Reset Stars", kTbResetStars, hints4);
810 fList->Add(hints4);
811 fList->Add(hints5);
812
813
814 TGHorizontalFrame *comp = new TGHorizontalFrame(grp2, 1, 1);
815 grp2->AddFrame(comp);
816 fList->Add(comp);
817
818 TGLayoutHints *hints3 = new TGLayoutHints(kLHintsLeft|kLHintsTop, 0, 20, 5, 0);
819 fList->Add(hints3);
820
821 TGLabel *l;
822
823 TGVerticalFrame *vframe = new TGVerticalFrame(comp, 1, 1);
824 comp->AddFrame(vframe, hints3);
825 fList->Add(vframe);
826
827 for (int i=0; i<MBending::GetNumPar(); i++)
828 AddCheckButton(vframe, fBending.GetName(i), i);
829
830 vframe = new TGVerticalFrame(comp, 1, 1);
831 comp->AddFrame(vframe, hints3);
832 fList->Add(vframe);
833
834 l = new TGLabel(vframe, "+000.0000");
835 l->SetTextJustify(kTextRight);
836 fList->Add(l);
837 fLabel.Add(l);
838
839 TGButton *but = (TGButton*)fList->FindWidget(0);
840
841 TGLayoutHints *h = new TGLayoutHints(kLHintsCenterY, 0, 0, but->GetHeight()-l->GetHeight());
842 fList->Add(h);
843
844 vframe->AddFrame(l,h);
845
846 for (int i=1; i<MBending::GetNumPar(); i++)
847 AddLabel(vframe, "+000.0000", h)->SetTextJustify(kTextRight);
848
849 vframe = new TGVerticalFrame(comp, 1, 1);
850 comp->AddFrame(vframe, hints3);
851 fList->Add(vframe);
852
853 for (int i=0; i<MBending::GetNumPar(); i++)
854 AddLabel(vframe, "\xb1 00.0000\xb0", h)->SetTextJustify(kTextRight);
855
856 vframe = new TGVerticalFrame(comp, 1, 1);
857 comp->AddFrame(vframe, hints3);
858 fList->Add(vframe);
859
860 for (int i=0; i<MBending::GetNumPar(); i++)
861 AddLabel(vframe, fBending.GetDescription(i), h);
862
863 l = new TGLabel(grp1, "0000000 Data Sets loaded.");
864 grp1->AddFrame(l, hints5);
865 fList->Add(l);
866 fLabel.Add(l);
867
868
869 ((TGCheckButton*)fList->FindWidget(0))->SetState(kButtonDown);
870 ((TGCheckButton*)fList->FindWidget(1))->SetState(kButtonDown);
871 ((TGCheckButton*)fList->FindWidget(3))->SetState(kButtonDown);
872 ((TGCheckButton*)fList->FindWidget(4))->SetState(kButtonDown);
873 ((TGCheckButton*)fList->FindWidget(5))->SetState(kButtonDown);
874 ((TGCheckButton*)fList->FindWidget(6))->SetState(kButtonDown);
875 ((TGCheckButton*)fList->FindWidget(7))->SetState(kButtonDown);
876
877 SetWindowName("TPoint Fitting Window");
878 SetIconName("TPoint++");
879
880 Layout();
881
882 MapSubwindows();
883 MapWindow();
884
885 DisplayBending();
886 DisplayData();
887 }
888 ClassDef(MFit, 0)
889};
890
891ClassImp(MFit);
892
893void gui()
894{
895 new MFit;
896}
Note: See TracBrowser for help on using the repository browser.