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