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