source: branches/Corsika7405Compatibility/mpointing/MPointingDevCalc.cc@ 18459

Last change on this file since 18459 was 14861, checked in by tbretz, 12 years ago
Updated all functions related to MAstro::GetDevAbs - they might however still need some check for signs.
File size: 27.7 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, 07/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2000-2007
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MPointingDevCalc
28//
29// Calculates the pointing deviation from the starguider information.
30//
31// There are some quality parameters to steer which are the valid
32// starguider data points:
33//
34// * MPointingDevCalc.NumMinStars: 8
35// Minimum number of identified stars required to accep the data
36//
37// * MPointingDevCalc.NsbLevel: 3.0
38// Minimum deviation (in rms) from the the mean allowed for the
39// measured NSB (noise in the ccd camera)
40//
41// * MPointingDevCalc.NsbMin: 30
42// - minimum NSB to be used in mean/rms calculation
43//
44// * MPointingDevCalc.NsbMax: 60
45// - maximum NSB to be used in mean/rms calculation
46//
47// * MPointingDevCalc.MaxAbsDev: 15
48// - Maximum absolute deviation which is consideres as valid (arcmin)
49//
50// Starguider data which doens't fullfill this requirements is ignored.
51// If the measures NSB==0 (file too old, starguider didn't write down
52// these values) the checks based on NSB and NumMinStar are skipped.
53//
54// The calculation of NSB mean and rms is reset for each file (ReInit)
55//
56// If your starguider data doesn't fullfill this requirement the latest
57// value which could be correctly calculated is used instead. If the time
58// for which no valid value can be calculated exceeds one minute
59// the return value is reset to 0/0. The maximum time allowed without
60// a valid value can be setup using:
61//
62// * MPointingDevCalc.MaxAge: 1
63// Maximum time before the starguider is reset to 0/0 in minutes
64//
65// Note, that the starguider itself is not well calibrated. Therefore
66// it is necessary to do a starguider calibration in our software.
67//
68// There are two options:
69//
70// * Simple starguider calibration using offsets in the camera plane
71//
72// The starguider is calibrated by taking its values (dZd/dAz)
73// adding them to the source position, calculating the source position
74// in the camera plane and adding the offsets. To switch off the
75// full starguider calibration do:
76//
77// * MPointingDevCalc.PointingModels:
78//
79// To set the offsets (in units of degree) use
80//
81// * MPointingDevCalc.Dx: -0.001
82// * MPointingDevCalc.Dy: -0.004
83//
84// * A starguider calibration using a pointing model calculated
85// from calibration data, so called TPoints
86//
87// Because the pointing model can change from time to time
88// you can give the run-number from which on a new pointing
89// model is valid. The run itself is included, e.g.:
90//
91// * MPointingDevCalc.PointingModels: 85240 89180
92// * MPointingDevCalc.FilePrefix: resources/starguider
93//
94// mean that for all runs<85240 the simple offset correction is used.
95// For runs >=85240 and <89180 the file resources/starguider0085240.txt
96// and for runs >=89180 the file resources/starguider0089180.txt is
97// used. To setup a default file for all runs before 85240 setup
98// a low number (eg. 0 or 1)
99//
100// In the case a pointing model is used additional offsets in
101// the x/y-camera plane (in units of deg) can be set using the DX
102// and DY parameters of the pointing model. The fDx and fDy data
103// members of this class are ignored. To overwrite the starguider
104// calibrated offset in either Az or Zd with a constant, you
105// can use the PX/PY directive in the pointing model. (To enable
106// the overwrite set the third column, the error, to a value
107// greater than zero)
108//
109//
110// At the PostProcessing step a table with statistics is print if the
111// debug level is greater or equal 3 (in most applications it is switched
112// on by -v3)
113//
114//
115// Pointing Models:
116// ----------------
117//
118// What we know so far about (maybe) important changes in cosy or to the
119// hardware:
120//
121// 25.01.2010: Changed scaling factors for TPoint camera (M1)
122//
123// 11.01.2010: New TPoint camera (M1)
124//
125// 18.03.2006: The camera holding has been repaired and the camera got
126// shifted a little bit.
127//
128// 16.04.2006: Around this date a roque lamp focussing has been done
129//
130// 25.04.2006: A missalignment introduced with the roque adjust has been
131// corrected
132//
133// The starguider pointing model for the time before 18.3.2006 and after
134// April 2006 (in fact there are no TPoints until 07/2006 to check it)
135// and for the period 07/2006 to (at least) 06/2007 are very similar.
136//
137// The pointing model for the time between 18.3.2006 and 04/2006 is
138// clearly different, mainly giving different Azimuth values between
139// Zenith and roughly ~25deg, and a slight offset on both axes.
140//
141// 10.5.2006: pos1 -= pos0 commented (What was the mentioned fix?)
142// 29.6.2006: repaired
143//
144// 23.3.2006: new starguider algorithm
145//
146// 17.3.2005: Fixed units of "nompos" in MDriveCom
147//
148//
149// New pointing models have been installed (if the pointing model
150// is different, than the previous one it might mean, that also
151// the starguider calibration is different.) The date only means
152// day-time (noon) at which the model has been changed.
153//
154// 29. Apr. 2004 ~25800
155// 5. Aug. 2004 ~32000
156// 19. Aug. 2004 ~33530
157// 7. Jun. 2005 ~57650
158// 8. Jun. 2005
159// 9. Jun. 2005 ~57860
160// 12. Sep. 2005 ~68338
161// 24. Nov. 2005 ~75562
162// 17. Oct. 2006 ~103130
163// 17. Jun. 2007 ~248193
164// 18. Oct. 2007 291104 // Correction for the offsets introduced by AMC
165// 14. Jan. 2008 328198 // Complete new pointing model
166// 11. Jun. 2008 1000534 (ca. 23:00) // Before new TPoints
167// 19. Jun. 2008 (ca. 15:00) // From new TPoints
168// 7. Mar. 2009 (ca. 14:00) // From new TPoints (0808-0902)
169// 14. May 2009 // M1/M2 after upgrade (from TPoints taken in the days before)
170// 17. Aug. 2009 New pointing models for both telescopes
171// 01. Feb. 2010 New pointing models for both telescopes
172// 26. Feb. 2010 New pointing models for both telescopes
173// 31. Mar. 2010 New pointing models for both telescopes
174// 29. Sep. 2010 New pointing model MAGIC I
175// 01. Dec. 2010 New pointing models for both telescopes
176// 06. Jan. 2011 New pointing models for both telescopes
177//
178// From 2.2.2006 beginnig of the night (21:05h, run >=81855) to 24.2.2006
179// beginning of the the night (20:34h, run<83722) the LEDs did not work.
180// In the time after this incident the shift crew often forgot to switch on
181// the LEDs at the beginning of the night!
182//
183// [2006-03-07 00:10:58] In the daytime, we raised the position of the
184// 9 o'clock LED by one screw hole to make it visible when the TPoint
185// Lid is closed. (< run 84980)
186//
187// Mirror refocussing around 23.Apr.2006, from the Runbook.
188//
189// 25./26.4.2006: Run 89180
190//
191// (Note: The year here is not a typo!)
192// [2007-04-25 23:50:39]
193// Markus is performing AMC focussing.
194//
195// Mirror refocussing around 4.Aug.2007, from the Runbook:
196//
197// [2007-08-04 04:46:47]
198// We finished with the focussing with Polaris. The images need to be
199// analysed and new LUTs generated.
200//
201// [2007-08-04 23:47:30]
202// Actually we see that the mispointing is always large; probably since
203// the LUT tables have not yet been adjusted to the new focussing.
204//
205// [2007-08-03 23:07:58]
206// Data taking stopped. Mirror focussing.
207//
208// [2007-08-05 00:09:16]
209// We take some pictures on stars nearby Cyg X3 with the sbig camera;
210// actually the spot doesn't look very nice... The pictures have been
211// saved with name Deneb- and Sadr- Polaris seems a bit better. Should we
212// have new LUT tables after the focussing?
213//
214// [2007-08-10 20:18:48]
215// Tonight we take first images of Polaris with a new LUT file generated
216// based on the recent focussing. The image will be analysed tomorrow and
217// than new LUTs will be generated. For tonight the focussing is still not
218// changed.
219//
220// [2007-08-14 20:57:59]
221// The weather is fine. There is a group of hobby astronomers at the
222// helicopter parking. Before data taking, we tried to check the new LUTs.
223// However, because of technical problems with the new LUTs we had to
224// postpone the measurements. We lost some 10 min of data taking because
225// of this.
226//
227// [2007-08-14 22:29:37] Run 267253
228// Before continuing the observation we perform a focussing test with the
229// new LUTs from recent Polaris focussing. Note: Data on Her X-1 was taken
230// with old focussing. We performed PSF measurements on Kornephorus (Zd
231// 31.77, Az 264,67) first with old LUTs and then with new LUTs. We took
232// T-points with both focussing. The first T-Point corresponds to the
233// previous focussing and the second with the improved. Please check both
234// T-points for eventual misspointing. We found big improvement of the PSF
235// with the new LUTs and will therefore continue from now on with the new
236// focussing. Having a first look at the SBIG pictures we see a slightly
237// misspointing of ~0.1 deg with the new LUTs.
238//
239// [2007-08-14 22:46:10]
240// Comparing the trigger rate with yesterday night we do not see an
241// improvement with the new focussing.
242//
243// [2007-10-27 email]
244// [...]
245// When removing the cover of the motor all 4 nuts which fix the gear of
246// the motor to the structure fell out. The motor was completely loose.
247// [...]
248// We checked also the second azimuth motor and .. the same situation.
249// [...]
250// Peter Sawallisch opened and fixed all 3 motors of MAGIC-I. We checked
251// all critical screws that came to our minds and fixed them as much as
252// possible with loctite and counter nuts. No damage whatsoever has been
253// found.
254//
255// [2009-08-05 Mail from Markus Garcz. about M1]
256// [...]
257// On the 23.07.2009 we tested the new LUTs using the star Etamin. The
258// PSF (sigma) improved from 11.0 mm to 9.5 mm and the total light content
259// in one PMT from 54% to 62%.
260// Furthermore I saw a small movement of the spot after the new focussing:
261// dX = 3 CCD pixel = 7 mm
262// dY = 2 CCD pixel = 4.7 mm
263// The new LUTs are now installed since the 23.07.2009 and are used as
264// default during the observation.
265// [...]
266//
267//
268// [2009-10-16 Email Adrian]
269//
270// I checked the AMC2 log files to check which LUTs had been used in which nights
271// until 09/09/12 the old version from March had always been used
272// 09/09/13-09/09/16 LUT_090913 (test) was used (by mistake)
273// 09/09/17-09/09/20 switched back to the version from March
274// since 09/09/21 using actual version LUT_090918
275// (PSF at ~10deg still rather poor, but I see strange
276// effects I do not yet understand and therefore cannot
277// correct)
278// Especially: nothing has changed in the AMC on 09/09/09,
279// except that around that date the SBIG camera was readjusted (but this
280// has no direct affect on the AMC)
281//
282//
283// [2010-02-03 M. Garcz.]
284//
285// New LUTs for M2
286//
287// [2010-06-14 M. Garcz.]
288//
289// New LUTs for M1
290//
291// [2010-07-02 M. Garcz.]
292//
293// New LUTs for M1
294//
295// [2010-10-24]
296//
297// Since about this date we see a higher spread in the TPoints. The reason we
298// found is that SA sends the coordinates of the orginal position and not the
299// target position to the AMC thus the mirrors are not well focussed.
300// This should mainly effect the quality of the TPoints but not the pointing
301// model. However, M2 shows a slight offset.
302//
303// [2010-11-11 M. Garcz.]
304//
305// Last night, TPoints for M2 taken without target
306//
307// [2010-11-18 M. Garcz.]
308//
309// New LUTs for M1
310//
311//
312// Others
313// ------
314//
315// The pointing of both(!) telescopes changed >2009/06/11 2:36h
316// In both cases the reason is unknown.
317//
318// The pointing has changed for M2 after 2009/09/09 noon
319// It recovered to the previous pointing at 2009/10/08 noon
320// Resons unknown (there were LUT changes but at different dates)
321//
322//
323// ToDo:
324// -----
325//
326// * Is 0/0 the best assumption if the starguider partly fails?
327//
328//
329// Input Container:
330// MRawRunHeader
331// MReportStarguider
332//
333// Output Container:
334// MPointingDev
335//
336/////////////////////////////////////////////////////////////////////////////
337#include "MPointingDevCalc.h"
338
339#include <stdlib.h> // atoi (Ubuntu 8.10)
340
341#include "MLog.h"
342#include "MLogManip.h"
343
344#include "MParList.h"
345
346#include "MAstro.h"
347#include "MPointing.h"
348#include "MPointingDev.h"
349#include "MRawRunHeader.h"
350#include "MReportStarguider.h"
351
352ClassImp(MPointingDevCalc);
353
354using namespace std;
355
356const TString MPointingDevCalc::fgFilePrefix="resources/starguider";
357
358// --------------------------------------------------------------------------
359//
360// Destructor. Call Clear() and delete fPointingModels if any.
361//
362MPointingDevCalc::~MPointingDevCalc()
363{
364 Clear();
365}
366
367// --------------------------------------------------------------------------
368//
369// Delete fPointing and set pointing to NULL
370//
371void MPointingDevCalc::Clear(Option_t *o)
372{
373 if (fPointing)
374 delete fPointing;
375
376 fPointing = NULL;
377}
378
379// --------------------------------------------------------------------------
380//
381// Sort the entries in fPoinitngModels
382//
383void MPointingDevCalc::SortPointingModels()
384{
385 const int n = fPointingModels.GetSize();
386
387 TArrayI idx(n);
388
389 TMath::Sort(n, fPointingModels.GetArray(), idx.GetArray(), kFALSE);
390
391 const TArrayI arr(fPointingModels);
392
393 for (int i=0; i<n; i++)
394 fPointingModels[i] = arr[idx[i]];
395}
396
397// --------------------------------------------------------------------------
398//
399// Set new pointing models
400//
401void MPointingDevCalc::SetPointingModels(const TString &models)
402{
403 fPointingModels.Set(0);
404
405 if (models.IsNull())
406 return;
407
408 TObjArray *arr = models.Tokenize(" ");
409
410 const int n = arr->GetEntries();
411 fPointingModels.Set(n);
412
413 for (int i=0; i<n; i++)
414 fPointingModels[i] = atoi((*arr)[i]->GetName());
415
416 delete arr;
417
418 SortPointingModels();
419}
420
421// --------------------------------------------------------------------------
422//
423// Return a string with the pointing models, seperated by a space.
424//
425TString MPointingDevCalc::GetPointingModels() const
426{
427 TString rc;
428 for (int i=0; i<fPointingModels.GetSize(); i++)
429 rc += Form ("%d ", fPointingModels[i]);
430
431 return rc;
432}
433
434// --------------------------------------------------------------------------
435//
436// Add a number to the pointing models
437//
438void MPointingDevCalc::AddPointingModel(UInt_t runnum)
439{
440 const int n = fPointingModels.GetSize();
441 for (int i=0; i<n; i++)
442 if ((UInt_t)fPointingModels[i]==runnum)
443 {
444 *fLog << warn << "WARNING - Pointing model " << runnum << " already in list... ignored." << endl;
445 return;
446 }
447
448 fPointingModels.Set(n+1);
449 fPointingModels[n] = runnum;
450
451 SortPointingModels();
452}
453
454// --------------------------------------------------------------------------
455//
456// Find the highest number in the array which is lower or equal num.
457//
458UInt_t MPointingDevCalc::FindPointingModel(UInt_t num)
459{
460 const int n = fPointingModels.GetSize();
461 if (n==0)
462 return (UInt_t)-1;
463
464 // Loop over all pointing models
465 for (int i=0; i<n; i++)
466 {
467 // The number stored in the array did not yet overtake the runnumber
468 if ((UInt_t)fPointingModels[i]<num)
469 continue;
470
471 // The first pointing model is later than this run: use a default
472 if (i==0)
473 return (UInt_t)-1;
474
475 // The last entry in the array is the right one.
476 return fPointingModels[i-1];
477 }
478
479 // Runnumber is after last entry of pointing models. Use the last one.
480 return fPointingModels[n-1];
481}
482
483// --------------------------------------------------------------------------
484//
485// Clear the pointing model. If run-number >= 87751 read the new
486// pointing model with fFilePrefix
487//
488Bool_t MPointingDevCalc::ReadPointingModel(const MRawRunHeader &run)
489{
490 const UInt_t num = FindPointingModel(run.GetRunNumber());
491
492 // No poinitng models are defined. Use simple dx/dy-calibration
493 if (num==(UInt_t)-1)
494 {
495 Clear();
496 return kTRUE;
497 }
498
499 // compile the name for the starguider files
500 // The file with the number 00000000 is the default file
501 TString fname = Form("%s%08d.txt", fFilePrefix.Data(), num);
502
503 if (!fPointing)
504 fPointing = new MPointing;
505
506 if (fname==fPointing->GetName())
507 {
508 *fLog << inf << fname << " already loaded." << endl;
509 return kTRUE;
510 }
511
512 return fPointing->Load(fname);
513}
514
515// --------------------------------------------------------------------------
516//
517// Check the file/run type from the run-header and if it is a data file
518// load starguider calibration.
519//
520Bool_t MPointingDevCalc::ReInit(MParList *plist)
521{
522 MRawRunHeader *run = (MRawRunHeader*)plist->FindObject("MRawRunHeader");
523 if (!run)
524 {
525 *fLog << err << "MRawRunHeader not found... aborting." << endl;
526 return kFALSE;
527 }
528
529 fNsbSum = 0;
530 fNsbSq = 0;
531 fNsbCount = 0;
532
533 fRunType = run->GetRunType();
534
535 switch (fRunType)
536 {
537 case MRawRunHeader::kRTData:
538 if (!fReport)
539 *fLog << warn << "MReportStarguider not found... skipped." << endl;
540 return ReadPointingModel(*run);
541
542 case MRawRunHeader::kRTMonteCarlo:
543 return kTRUE;
544
545 case MRawRunHeader::kRTPedestal:
546 *fLog << err << "Cannot work in a pedestal Run!... aborting." << endl;
547 return kFALSE;
548
549 case MRawRunHeader::kRTCalibration:
550 *fLog << err << "Cannot work in a calibration Run!... aborting." << endl;
551 return kFALSE;
552
553 default:
554 *fLog << err << "Run Type " << fRunType << " unknown!... aborting." << endl;
555 return kFALSE;
556 }
557
558 return kTRUE;
559}
560
561// --------------------------------------------------------------------------
562//
563// Search for 'MPointingPos'. Create if not found.
564//
565Int_t MPointingDevCalc::PreProcess(MParList *plist)
566{
567 fDeviation = (MPointingDev*)plist->FindCreateObj("MPointingDev");
568 fReport = (MReportStarguider*)plist->FindObject("MReportStarguider");
569
570 // We use kRTNone here as a placeholder for data runs.
571 fRunType = MRawRunHeader::kRTNone;
572 fLastMjd = -1;
573
574 fSkip.Reset();
575
576 // In cases in which ReInit isn't called right in time (e.g. if
577 // the first starguider event comes before the first data event)
578 fNsbSum = 0;
579 fNsbSq = 0;
580 fNsbCount = 0;
581
582 return fDeviation ? kTRUE : kFALSE;
583}
584
585// --------------------------------------------------------------------------
586//
587// Increase fSkip[i] by one. If the data in fDeviation is outdated (older
588// than fMaxAge) and the current report should be skipped reset DevZdAz and
589// DevXY and fSkip[6] is increased by one.
590//
591void MPointingDevCalc::Skip(Int_t i)
592{
593 fSkip[i]++;
594
595 const Double_t diff = (fReport->GetMjd()-fLastMjd)*1440; // [min] 1440=24*60
596 if (diff<fMaxAge && fLastMjd>0)
597 return;
598
599 fDeviation->SetDevZdAz(0, 0);
600 fDeviation->SetDevXY(0, 0);
601 fSkip[6]++;
602}
603
604// --------------------------------------------------------------------------
605//
606// Do a full starguider calibration using a pointing model for the starguider.
607//
608void MPointingDevCalc::DoCalibration(Double_t devzd, Double_t devaz) const
609{
610 if (!fPointing)
611 {
612 // Do a simple starguider calibration using a simple offset in x and y
613 fDeviation->SetDevZdAz(devzd, devaz);
614
615 // Linear starguider calibration taken from April/May data
616 // For calibration add MDriveReport::GetErrorZd/Az !
617 fDeviation->SetDevXY(fDx, fDy); // 1arcmin ~ 5mm
618
619 return;
620 }
621
622 // Get the nominal position the star is at the sky
623 // Unit: deg
624 ZdAz nom(fReport->GetNominalZd(), fReport->GetNominalAz());
625 nom *= TMath::DegToRad();
626
627 // Get the mispointing measured by the telescope. It is
628 // calculate as follows:
629 //
630 // The mispointing measured by the starguider:
631 // ZdAz mis(devzd, devaz);
632 // mis *= TMath::DegToRad();
633
634 // The pointing model is the conversion from the real pointing
635 // position of the telescope into the pointing position measured
636 // by the starguider.
637 //
638 // To keep as close to the fitted model we use the forward correction.
639
640 // Position at which the starguider camera is pointing in real:
641 // pointing position = nominal position - dev
642 //
643 // The position measured as the starguider's pointing position
644 ZdAz pos(nom); // cpos = sao - dev
645 pos -= ZdAz(devzd, devaz)*TMath::DegToRad();
646
647 // Now we convert the starguider's pointing position into the
648 // telescope pointing position (the pointing model converts
649 // the telescope pointing position into the starguider pointing
650 // position)
651 ZdAz point = fPointing->CorrectBack(pos); //FWD!!!
652
653 // MSrcPosCalc uses the following condition to calculate the
654 // source position in the camera:
655 // real pointing pos = nominal pointing pos - dev
656 // --> dev = nominal - real
657 // Therefor we calculate dev as follows:
658 ZdAz dev(nom);
659 dev -= point;
660 dev *= TMath::RadToDeg();
661
662 /*
663 // We chose the other way. It is less accurate because is is the
664 // other was than the poinitng model was fittet, but it is more
665 // accurate because the nominal (i.e. real) pointing position
666 // is less accurately known than the position returned by the
667 // starguider.
668 //
669 // Calculate the deviation which would be measured by the starguider
670 // if applied to a perfectly pointing telescope.
671 ZdAz dev = fPointing->Correct(nom);
672 dev -= nom;
673
674 // Now add these offsets and the starguider measured offsets to
675 // the real pointing deviation of the telescope (note, that
676 // signs here are just conventions)
677 dev += ZdAz(devzd, devaz)*TMath::DegToRad(); // --> nom-mis
678 dev *= TMath::RadToDeg();
679 */
680
681 // Check if the starguider pointing model requests overwriting
682 // of the values with constants (e.g. 0)
683 devaz = fPointing->IsPxValid() ? fPointing->GetPx() : dev.Az();
684 devzd = fPointing->IsPyValid() ? fPointing->GetPy() : dev.Zd();
685
686 fDeviation->SetDevZdAz(devzd, devaz);
687 fDeviation->SetDevXY(fPointing->GetDxy());
688}
689
690Int_t MPointingDevCalc::ProcessStarguiderReport()
691{
692 Double_t devzd = fReport->GetDevZd(); // [arcmin]
693 Double_t devaz = fReport->GetDevAz(); // [arcmin]
694 if (devzd==0 && devaz==0)
695 {
696 Skip(1);
697 return kTRUE;
698 }
699
700 if (!fReport->IsMonitoring())
701 {
702 Skip(2);
703 return kTRUE;
704 }
705
706 devzd /= 60; // Convert from arcmin to deg
707 devaz /= 60; // Convert from arcmin to deg
708
709 const Double_t nsb = fReport->GetSkyBrightness();
710 if (nsb>0)
711 {
712 if (nsb>fNsbMin && nsb<fNsbMax)
713 {
714 fNsbSum += nsb;
715 fNsbSq += nsb*nsb;
716 fNsbCount++;
717 }
718
719 if (fNsbCount>0)
720 {
721 const Double_t sum = fNsbSum/fNsbCount;
722 const Double_t sq = fNsbSq /fNsbCount;
723
724 const Double_t rms = fNsbLevel*TMath::Sqrt(sq - sum*sum);
725
726 if (nsb<sum-rms || nsb>sum+rms)
727 {
728 Skip(3);
729 return kTRUE;
730 }
731 }
732
733 if (fReport->GetNumIdentifiedStars()<fNumMinStars)
734 {
735 Skip(4);
736 return kTRUE;
737 }
738 }
739
740 // >= 87751 (31.3.06 12:00)
741
742 // Calculate absolute deviation
743 const Double_t dev = MAstro::GetDevAbs(fReport->GetNominalZd(),
744 fReport->GetNominalZd()-fReport->GetDevZd()/60,
745 devaz);
746
747 // Sanity check... larger deviation are strange and ignored
748 if (dev*60>fMaxAbsDev)
749 {
750 Skip(5);
751 return kTRUE;
752 }
753
754 DoCalibration(devzd, devaz);
755
756 fSkip[0]++;
757 fLastMjd = fReport->GetMjd();
758
759 return kTRUE;
760}
761
762// --------------------------------------------------------------------------
763//
764// See class description.
765//
766Int_t MPointingDevCalc::Process()
767{
768 switch (fRunType)
769 {
770 case MRawRunHeader::kRTNone:
771 case MRawRunHeader::kRTData:
772 return fReport ? ProcessStarguiderReport() : kTRUE;
773
774 case MRawRunHeader::kRTMonteCarlo:
775 fSkip[0]++;
776 fDeviation->SetDevZdAz(0, 0);
777 fDeviation->SetDevXY(0, 0);
778 return kTRUE;
779 }
780 return kTRUE;
781}
782
783// --------------------------------------------------------------------------
784//
785// Print execution statistics
786//
787Int_t MPointingDevCalc::PostProcess()
788{
789 if (GetNumExecutions()==0)
790 return kTRUE;
791
792 *fLog << inf << endl;
793 *fLog << GetDescriptor() << " execution statistics:" << endl;
794 PrintSkipped(fSkip[1], "Starguider deviation not set, is exactly 0/0");
795 PrintSkipped(fSkip[2], "Starguider was not monitoring (eg. LEDs off)");
796 PrintSkipped(fSkip[3], Form("NSB out of %.1f sigma range", fNsbLevel));
797 PrintSkipped(fSkip[4], Form("Number of identified stars < %d", fNumMinStars));
798 PrintSkipped(fSkip[5], Form("Absolute deviation > %.1farcmin", fMaxAbsDev));
799 PrintSkipped(fSkip[6], Form("Events set to 0 because older than %.1fmin", fMaxAge));
800 *fLog << " " << (int)fSkip[0] << " (" << Form("%5.1f", 100.*fSkip[0]/GetNumExecutions()) << "%) Evts survived calculation!" << endl;
801 *fLog << endl;
802
803 return kTRUE;
804}
805
806// --------------------------------------------------------------------------
807//
808// MPointingDevCalc.NumMinStars: 8
809// MPointingDevCalc.NsbLevel: 3.0
810// MPointingDevCalc.NsbMin: 30
811// MPointingDevCalc.NsbMax: 60
812// MPointingDevCalc.MaxAbsDev: 15
813// MPointingDevCalc.MaxAge: 1.0
814// MPointingDevCalc.Dx: -0.001
815// MPointingDevCalc.Dy: -0.004
816//
817// For a detailed description see the class reference.
818//
819Int_t MPointingDevCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
820{
821 Bool_t rc = kFALSE;
822 if (IsEnvDefined(env, prefix, "NumMinStars", print))
823 {
824 SetNumMinStars(GetEnvValue(env, prefix, "NumMinStars", (Int_t)fNumMinStars));
825 rc = kTRUE;
826 }
827 if (IsEnvDefined(env, prefix, "NsbLevel", print))
828 {
829 SetNsbLevel(GetEnvValue(env, prefix, "NsbLevel", fNsbLevel));
830 rc = kTRUE;
831 }
832 if (IsEnvDefined(env, prefix, "NsbMin", print))
833 {
834 SetNsbMin(GetEnvValue(env, prefix, "NsbMin", fNsbMin));
835 rc = kTRUE;
836 }
837 if (IsEnvDefined(env, prefix, "NsbMax", print))
838 {
839 SetNsbMax(GetEnvValue(env, prefix, "NsbMax", fNsbMax));
840 rc = kTRUE;
841 }
842 if (IsEnvDefined(env, prefix, "MaxAbsDev", print))
843 {
844 SetMaxAbsDev(GetEnvValue(env, prefix, "MaxAbsDev", fMaxAbsDev));
845 rc = kTRUE;
846 }
847 if (IsEnvDefined(env, prefix, "Dx", print))
848 {
849 fDx = GetEnvValue(env, prefix, "Dx", fDx);
850 rc = kTRUE;
851 }
852 if (IsEnvDefined(env, prefix, "Dy", print))
853 {
854 fDy = GetEnvValue(env, prefix, "Dy", fDy);
855 rc = kTRUE;
856 }
857 if (IsEnvDefined(env, prefix, "MaxAge", print))
858 {
859 fMaxAge = GetEnvValue(env, prefix, "MaxAge", fMaxAge);
860 rc = kTRUE;
861 }
862 if (IsEnvDefined(env, prefix, "FilePrefix", print))
863 {
864 fFilePrefix = GetEnvValue(env, prefix, "FilePrefix", fFilePrefix);
865 rc = kTRUE;
866 }
867 if (IsEnvDefined(env, prefix, "PointingModels", print))
868 {
869 SetPointingModels(GetEnvValue(env, prefix, "PointingModels", GetPointingModels()));
870 rc = kTRUE;
871 }
872
873 return rc;
874}
Note: See TracBrowser for help on using the repository browser.