Commit d495224164d0efdeabe1e37f384267b83c7ee714
1 parent
eb128f34
initial version
Showing
5 changed files
with
822 additions
and
0 deletions
src/ctskymap/Makefile.am
0 → 100644
| 1 | +## Process this file with automake to produce Makefile.in | ||
| 2 | + | ||
| 3 | +# Program name | ||
| 4 | +bin_PROGRAMS = ctskymap | ||
| 5 | + | ||
| 6 | +# Sources | ||
| 7 | +sources = ctskymap.cpp | ||
| 8 | + | ||
| 9 | +# Libtool library | ||
| 10 | +noinst_LTLIBRARIES = libctskymap.la | ||
| 11 | +libctskymap_la_SOURCES = $(sources) | ||
| 12 | + | ||
| 13 | +# Program definition | ||
| 14 | +ctskymap_SOURCES = main.cpp | ||
| 15 | +ctskymap_LDFLAGS = $(LIB_LINK_OPTION) @LDFLAGS@ | ||
| 16 | +ctskymap_CPPFLAGS = @CPPFLAGS@ | ||
| 17 | +ctskymap_LDADD = libctskymap.la | ||
| 18 | + | ||
| 19 | +# Parameter files | ||
| 20 | +pfilesdir = $(prefix)/syspfiles | ||
| 21 | +pfiles_DATA = ctskymap.par |
src/ctskymap/ctskymap.cpp
0 → 100644
| 1 | +/*************************************************************************** | ||
| 2 | + * ctskymap - CTA sky mapping tool * | ||
| 3 | + * ----------------------------------------------------------------------- * | ||
| 4 | + * copyright (C) 2010-2011 by Jurgen Knodlseder * | ||
| 5 | + * ----------------------------------------------------------------------- * | ||
| 6 | + * * | ||
| 7 | + * This program is free software: you can redistribute it and/or modify * | ||
| 8 | + * it under the terms of the GNU General Public License as published by * | ||
| 9 | + * the Free Software Foundation, either version 3 of the License, or * | ||
| 10 | + * (at your option) any later version. * | ||
| 11 | + * * | ||
| 12 | + * This program is distributed in the hope that it will be useful, * | ||
| 13 | + * but WITHOUT ANY WARRANTY; without even the implied warranty of * | ||
| 14 | + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * | ||
| 15 | + * GNU General Public License for more details. * | ||
| 16 | + * * | ||
| 17 | + * You should have received a copy of the GNU General Public License * | ||
| 18 | + * along with this program. If not, see <http://www.gnu.org/licenses/>. * | ||
| 19 | + * * | ||
| 20 | + ***************************************************************************/ | ||
| 21 | +/** | ||
| 22 | + * @file ctskymap.cpp | ||
| 23 | + * @brief CTA sky mapping tool implementation | ||
| 24 | + * @author J. Knodlseder | ||
| 25 | + */ | ||
| 26 | + | ||
| 27 | +/* __ Includes ___________________________________________________________ */ | ||
| 28 | +#ifdef HAVE_CONFIG_H | ||
| 29 | +#include <config.h> | ||
| 30 | +#endif | ||
| 31 | +#include <cstdio> | ||
| 32 | +#include "ctskymap.hpp" | ||
| 33 | +#include "GTools.hpp" | ||
| 34 | + | ||
| 35 | +/* __ Method name definitions ____________________________________________ */ | ||
| 36 | +#define G_BIN_EVENTS "ctskymap::bin_events(GCTAObservation*)" | ||
| 37 | + | ||
| 38 | +/* __ Debug definitions __________________________________________________ */ | ||
| 39 | + | ||
| 40 | +/* __ Coding definitions _________________________________________________ */ | ||
| 41 | + | ||
| 42 | + | ||
| 43 | +/*========================================================================== | ||
| 44 | + = = | ||
| 45 | + = Constructors/destructors = | ||
| 46 | + = = | ||
| 47 | + ==========================================================================*/ | ||
| 48 | + | ||
| 49 | +/***********************************************************************//** | ||
| 50 | + * @brief Void constructor | ||
| 51 | + ***************************************************************************/ | ||
| 52 | +ctskymap::ctskymap(void) : GApplication(CTSKYMAP_NAME, CTSKYMAP_VERSION) | ||
| 53 | +{ | ||
| 54 | + // Initialise members | ||
| 55 | + init_members(); | ||
| 56 | + | ||
| 57 | + // Write header into logger | ||
| 58 | + log_header(); | ||
| 59 | + | ||
| 60 | + // Return | ||
| 61 | + return; | ||
| 62 | +} | ||
| 63 | + | ||
| 64 | + | ||
| 65 | +/***********************************************************************//** | ||
| 66 | + * @brief Observations constructor | ||
| 67 | + * | ||
| 68 | + * This method creates an instance of the class by copying an existing | ||
| 69 | + * observations container. | ||
| 70 | + ***************************************************************************/ | ||
| 71 | +ctskymap::ctskymap(GObservations obs) : GApplication(CTSKYMAP_NAME, CTSKYMAP_VERSION) | ||
| 72 | +{ | ||
| 73 | + // Initialise members | ||
| 74 | + init_members(); | ||
| 75 | + | ||
| 76 | + // Set observations | ||
| 77 | + m_obs = obs; | ||
| 78 | + | ||
| 79 | + // Write header into logger | ||
| 80 | + log_header(); | ||
| 81 | + | ||
| 82 | + // Return | ||
| 83 | + return; | ||
| 84 | +} | ||
| 85 | + | ||
| 86 | + | ||
| 87 | + | ||
| 88 | +/***********************************************************************//** | ||
| 89 | + * @brief Command line constructor | ||
| 90 | + * | ||
| 91 | + * @param[in] argc Number of arguments in command line. | ||
| 92 | + * @param[in] argv Array of command line arguments. | ||
| 93 | + ***************************************************************************/ | ||
| 94 | +ctskymap::ctskymap(int argc, char *argv[]) : | ||
| 95 | + GApplication(CTSKYMAP_NAME, CTSKYMAP_VERSION, argc, argv) | ||
| 96 | +{ | ||
| 97 | + // Initialise members | ||
| 98 | + init_members(); | ||
| 99 | + | ||
| 100 | + // Write header into logger | ||
| 101 | + log_header(); | ||
| 102 | + | ||
| 103 | + // Return | ||
| 104 | + return; | ||
| 105 | +} | ||
| 106 | + | ||
| 107 | + | ||
| 108 | +/***********************************************************************//** | ||
| 109 | + * @brief Copy constructor | ||
| 110 | + * | ||
| 111 | + * @param[in] app Application. | ||
| 112 | + ***************************************************************************/ | ||
| 113 | +ctskymap::ctskymap(const ctskymap& app) : GApplication(app) | ||
| 114 | +{ | ||
| 115 | + // Initialise members | ||
| 116 | + init_members(); | ||
| 117 | + | ||
| 118 | + // Copy members | ||
| 119 | + copy_members(app); | ||
| 120 | + | ||
| 121 | + // Return | ||
| 122 | + return; | ||
| 123 | +} | ||
| 124 | + | ||
| 125 | + | ||
| 126 | +/***********************************************************************//** | ||
| 127 | + * @brief Destructor | ||
| 128 | + ***************************************************************************/ | ||
| 129 | +ctskymap::~ctskymap(void) | ||
| 130 | +{ | ||
| 131 | + // Free members | ||
| 132 | + free_members(); | ||
| 133 | + | ||
| 134 | + // Return | ||
| 135 | + return; | ||
| 136 | +} | ||
| 137 | + | ||
| 138 | + | ||
| 139 | +/*========================================================================== | ||
| 140 | + = = | ||
| 141 | + = Operators = | ||
| 142 | + = = | ||
| 143 | + ==========================================================================*/ | ||
| 144 | + | ||
| 145 | +/***********************************************************************//** | ||
| 146 | + * @brief Assignment operator | ||
| 147 | + * | ||
| 148 | + * @param[in] app Application. | ||
| 149 | + ***************************************************************************/ | ||
| 150 | +ctskymap& ctskymap::operator= (const ctskymap& app) | ||
| 151 | +{ | ||
| 152 | + // Execute only if object is not identical | ||
| 153 | + if (this != &app) { | ||
| 154 | + | ||
| 155 | + // Copy base class members | ||
| 156 | + this->GApplication::operator=(app); | ||
| 157 | + | ||
| 158 | + // Free members | ||
| 159 | + free_members(); | ||
| 160 | + | ||
| 161 | + // Initialise members | ||
| 162 | + init_members(); | ||
| 163 | + | ||
| 164 | + // Copy members | ||
| 165 | + copy_members(app); | ||
| 166 | + | ||
| 167 | + } // endif: object was not identical | ||
| 168 | + | ||
| 169 | + // Return this object | ||
| 170 | + return *this; | ||
| 171 | +} | ||
| 172 | + | ||
| 173 | + | ||
| 174 | +/*========================================================================== | ||
| 175 | + = = | ||
| 176 | + = Public methods = | ||
| 177 | + = = | ||
| 178 | + ==========================================================================*/ | ||
| 179 | + | ||
| 180 | +/***********************************************************************//** | ||
| 181 | + * @brief Clear instance | ||
| 182 | + ***************************************************************************/ | ||
| 183 | +void ctskymap::clear(void) | ||
| 184 | +{ | ||
| 185 | + // Free members | ||
| 186 | + free_members(); | ||
| 187 | + this->GApplication::free_members(); | ||
| 188 | + | ||
| 189 | + // Initialise members | ||
| 190 | + this->GApplication::init_members(); | ||
| 191 | + init_members(); | ||
| 192 | + | ||
| 193 | + // Return | ||
| 194 | + return; | ||
| 195 | +} | ||
| 196 | + | ||
| 197 | + | ||
| 198 | +/***********************************************************************//** | ||
| 199 | + * @brief Execute application | ||
| 200 | + * | ||
| 201 | + * This method bins the events data into a counts maps and saves the counts | ||
| 202 | + * map into a FITS file. | ||
| 203 | + ***************************************************************************/ | ||
| 204 | +void ctskymap::execute(void) | ||
| 205 | +{ | ||
| 206 | + // Read ahead output filename so that it gets dumped correctly in the | ||
| 207 | + // parameters log | ||
| 208 | + m_outfile = (*this)["outfile"].filename(); | ||
| 209 | + | ||
| 210 | + // Bin the event data | ||
| 211 | + run(); | ||
| 212 | + | ||
| 213 | + // Save the counts map into FITS file | ||
| 214 | + save(); | ||
| 215 | + | ||
| 216 | + // Return | ||
| 217 | + return; | ||
| 218 | +} | ||
| 219 | + | ||
| 220 | + | ||
| 221 | +/***********************************************************************//** | ||
| 222 | + * @brief Bin the event data | ||
| 223 | + ***************************************************************************/ | ||
| 224 | +void ctskymap::run(void) | ||
| 225 | +{ | ||
| 226 | + // Switch screen logging on in debug mode | ||
| 227 | + if (logDebug()) | ||
| 228 | + log.cout(true); | ||
| 229 | + | ||
| 230 | + // Get parameters | ||
| 231 | + get_parameters(); | ||
| 232 | + | ||
| 233 | + // Write parameters into logger | ||
| 234 | + if (logTerse()) { | ||
| 235 | + log_parameters(); | ||
| 236 | + log << std::endl; | ||
| 237 | + } | ||
| 238 | + | ||
| 239 | + // Write observation(s) into logger | ||
| 240 | + if (logTerse()) { | ||
| 241 | + log << std::endl; | ||
| 242 | + if (m_obs.size() > 1) | ||
| 243 | + log.header1("Observations"); | ||
| 244 | + else | ||
| 245 | + log.header1("Observation"); | ||
| 246 | + log << m_obs << std::endl; | ||
| 247 | + } | ||
| 248 | + | ||
| 249 | + // Write header | ||
| 250 | + if (logTerse()) { | ||
| 251 | + log << std::endl; | ||
| 252 | + if (m_obs.size() > 1) | ||
| 253 | + log.header1("Bin observations"); | ||
| 254 | + else | ||
| 255 | + log.header1("Bin observation"); | ||
| 256 | + } | ||
| 257 | + | ||
| 258 | + // Loop over all observation in the container | ||
| 259 | + for (int i = 0; i < m_obs.size(); ++i) { | ||
| 260 | + | ||
| 261 | + // Get CTA observation | ||
| 262 | + GCTAObservation* obs = dynamic_cast<GCTAObservation*>(&m_obs[i]); | ||
| 263 | + | ||
| 264 | + // Continue only if observation is a CTA observation | ||
| 265 | + if (obs != NULL) { | ||
| 266 | + | ||
| 267 | + // Write header for observation | ||
| 268 | + if (logTerse()) { | ||
| 269 | + if (obs->name().length() > 1) | ||
| 270 | + log.header3("Observation "+obs->name()); | ||
| 271 | + else | ||
| 272 | + log.header3("Observation"); | ||
| 273 | + } | ||
| 274 | + | ||
| 275 | + // Map events into sky map | ||
| 276 | + map_events(obs); | ||
| 277 | + | ||
| 278 | + } // endif: CTA observation found | ||
| 279 | + | ||
| 280 | + } // endfor: looped over observations | ||
| 281 | + | ||
| 282 | + // Write observation(s) into logger | ||
| 283 | + if (logTerse()) { | ||
| 284 | + log << std::endl; | ||
| 285 | + if (m_obs.size() > 1) | ||
| 286 | + log.header1("Binned observations"); | ||
| 287 | + else | ||
| 288 | + log.header1("Binned observation"); | ||
| 289 | + log << m_obs << std::endl; | ||
| 290 | + } | ||
| 291 | + | ||
| 292 | + // Return | ||
| 293 | + return; | ||
| 294 | +} | ||
| 295 | + | ||
| 296 | + | ||
| 297 | +/***********************************************************************//** | ||
| 298 | + * @brief Save observation | ||
| 299 | + * | ||
| 300 | + * This method saves the counts map(s) into (a) FITS file(s). | ||
| 301 | + ***************************************************************************/ | ||
| 302 | +void ctskymap::save(void) | ||
| 303 | +{ | ||
| 304 | +/* | ||
| 305 | + // Write header | ||
| 306 | + if (logTerse()) { | ||
| 307 | + log << std::endl; | ||
| 308 | + if (m_obs.size() > 1) | ||
| 309 | + log.header1("Save observations"); | ||
| 310 | + else | ||
| 311 | + log.header1("Save observation"); | ||
| 312 | + } | ||
| 313 | + | ||
| 314 | + // Get output filename | ||
| 315 | + m_outfile = (*this)["outfile"].filename(); | ||
| 316 | + | ||
| 317 | + // Loop over all observations in the container | ||
| 318 | + int file_num = 0; | ||
| 319 | + for (int i = 0; i < m_obs.size(); ++i) { | ||
| 320 | + | ||
| 321 | + // Get CTA observation | ||
| 322 | + GCTAObservation* obs = dynamic_cast<GCTAObservation*>(&m_obs[i]); | ||
| 323 | + | ||
| 324 | + // Save only if observation is a CTA observation | ||
| 325 | + if (obs != NULL) { | ||
| 326 | + | ||
| 327 | + // Set filename. If more than one file will be created an | ||
| 328 | + // index "_xxx" will be appended. | ||
| 329 | + std::string filename = m_outfile; | ||
| 330 | + if (file_num > 0) | ||
| 331 | + filename += "_"+str(file_num); | ||
| 332 | + | ||
| 333 | + // Save file | ||
| 334 | + obs->save(filename, clobber()); | ||
| 335 | + | ||
| 336 | + // Increment file number | ||
| 337 | + file_num++; | ||
| 338 | + | ||
| 339 | + } // endif: observation was a CTA observation | ||
| 340 | + | ||
| 341 | + } // endfor: looped over files | ||
| 342 | +*/ | ||
| 343 | + // Return | ||
| 344 | + return; | ||
| 345 | +} | ||
| 346 | + | ||
| 347 | + | ||
| 348 | +/***********************************************************************//** | ||
| 349 | + * @brief Get application parameters | ||
| 350 | + * | ||
| 351 | + * Get all task parameters from parameter file or (if required) by querying | ||
| 352 | + * the user. Most parameters are only required if no observation exists so | ||
| 353 | + * far in the observation container. In this case, a single CTA observation | ||
| 354 | + * will be added to the container, using the definition provided in the | ||
| 355 | + * parameter file. | ||
| 356 | + ***************************************************************************/ | ||
| 357 | +void ctskymap::get_parameters(void) | ||
| 358 | +{ | ||
| 359 | + // If there are no observations in container then add a single CTA | ||
| 360 | + // observation using the parameters from the parameter file | ||
| 361 | + if (m_obs.size() == 0) { | ||
| 362 | + | ||
| 363 | + // Get name of CTA events file | ||
| 364 | + m_evfile = (*this)["evfile"].filename(); | ||
| 365 | + | ||
| 366 | + // Load unbinned CTA observation | ||
| 367 | + GCTAObservation obs; | ||
| 368 | + obs.load_unbinned(m_evfile); | ||
| 369 | + | ||
| 370 | + // Append CTA observation to container | ||
| 371 | + m_obs.append(obs); | ||
| 372 | + | ||
| 373 | + // Use the xref and yref parameters for binning (otherwise the | ||
| 374 | + // pointing direction(s) is/are used) | ||
| 375 | + m_xref = (*this)["xref"].real(); | ||
| 376 | + m_yref = (*this)["yref"].real(); | ||
| 377 | + | ||
| 378 | + } // endif: there was no observation in the container | ||
| 379 | + | ||
| 380 | + // Get remaining parameters | ||
| 381 | + m_emin = (*this)["emin"].real(); | ||
| 382 | + m_emax = (*this)["emax"].real(); | ||
| 383 | + m_proj = (*this)["proj"].string(); | ||
| 384 | + m_coordsys = (*this)["coordsys"].string(); | ||
| 385 | + m_binsz = (*this)["binsz"].real(); | ||
| 386 | + m_nxpix = (*this)["nxpix"].integer(); | ||
| 387 | + m_nypix = (*this)["nypix"].integer(); | ||
| 388 | + | ||
| 389 | + // Return | ||
| 390 | + return; | ||
| 391 | +} | ||
| 392 | + | ||
| 393 | + | ||
| 394 | +/***********************************************************************//** | ||
| 395 | + * @brief Map events into a sky map | ||
| 396 | + * | ||
| 397 | + * @param[in] obs CTA observation. | ||
| 398 | + * | ||
| 399 | + * @exception GException::no_list | ||
| 400 | + * No event list found in observation. | ||
| 401 | + * @exception GCTAException::no_pointing | ||
| 402 | + * No valid CTA pointing found. | ||
| 403 | + * | ||
| 404 | + * This method maps the events found in a CTA events list into a sky map. | ||
| 405 | + ***************************************************************************/ | ||
| 406 | +void ctskymap::map_events(GCTAObservation* obs) | ||
| 407 | +{ | ||
| 408 | + // Continue only if observation pointer is valid | ||
| 409 | + if (obs != NULL) { | ||
| 410 | + | ||
| 411 | + // Make sure that the observation holds a CTA event list. If this | ||
| 412 | + // is not the case then throw an exception. | ||
| 413 | + if (dynamic_cast<const GCTAEventList*>(obs->events()) == NULL) { | ||
| 414 | + throw GException::no_list(G_BIN_EVENTS); | ||
| 415 | + } | ||
| 416 | + | ||
| 417 | + // Setup energy range covered by data | ||
| 418 | + GEnergy emin; | ||
| 419 | + GEnergy emax; | ||
| 420 | + GEbounds ebds; | ||
| 421 | + emin.TeV(m_emin); | ||
| 422 | + emax.TeV(m_emax); | ||
| 423 | + | ||
| 424 | + // Get map centre | ||
| 425 | + double xref; | ||
| 426 | + double yref; | ||
| 427 | + if (m_xref != 9999.0 && m_yref != 9999.0) { | ||
| 428 | + xref = m_xref; | ||
| 429 | + yref = m_yref; | ||
| 430 | + } | ||
| 431 | + else { | ||
| 432 | + | ||
| 433 | + // Get pointer on CTA pointing | ||
| 434 | + const GTime srcTime; // Dummy time | ||
| 435 | + const GCTAPointing *pnt = obs->pointing(srcTime); | ||
| 436 | + if (pnt == NULL) | ||
| 437 | + throw GCTAException::no_pointing(G_BIN_EVENTS); | ||
| 438 | + | ||
| 439 | + // Set reference point to pointing | ||
| 440 | + if (toupper(m_coordsys) == "GAL") { | ||
| 441 | + xref = pnt->dir().l_deg(); | ||
| 442 | + yref = pnt->dir().b_deg(); | ||
| 443 | + } | ||
| 444 | + else { | ||
| 445 | + xref = pnt->dir().ra_deg(); | ||
| 446 | + yref = pnt->dir().dec_deg(); | ||
| 447 | + } | ||
| 448 | + | ||
| 449 | + } // endelse: map centre set to pointing | ||
| 450 | + | ||
| 451 | + // Create skymap | ||
| 452 | + GSkymap map = GSkymap(m_proj, m_coordsys, | ||
| 453 | + xref, yref, m_binsz, m_binsz, | ||
| 454 | + m_nxpix, m_nypix, 1); | ||
| 455 | + | ||
| 456 | + // Initialise binning statistics | ||
| 457 | + int num_outside_map = 0; | ||
| 458 | + int num_outside_erange = 0; | ||
| 459 | + int num_in_map = 0; | ||
| 460 | + | ||
| 461 | + // Fill sky map | ||
| 462 | + GCTAEventList* events = (GCTAEventList*)obs->events(); | ||
| 463 | + for (GCTAEventList::iterator event = events->begin(); event != events->end(); ++event) { | ||
| 464 | + | ||
| 465 | + // Skip if energy is out of range | ||
| 466 | + if (event->energy() < emin || event->energy() > emax) { | ||
| 467 | + num_outside_erange++; | ||
| 468 | + continue; | ||
| 469 | + } | ||
| 470 | + | ||
| 471 | + // Determine sky pixel | ||
| 472 | + GCTAInstDir* inst = (GCTAInstDir*)&(event->dir()); | ||
| 473 | + GSkyDir dir = inst->skydir(); | ||
| 474 | + GSkyPixel pixel = map.dir2xy(dir); | ||
| 475 | + | ||
| 476 | + // Skip if pixel is out of range | ||
| 477 | + if (pixel.x() < -0.5 || pixel.x() > (m_nxpix-0.5) || | ||
| 478 | + pixel.y() < -0.5 || pixel.y() > (m_nypix-0.5)) { | ||
| 479 | + num_outside_map++; | ||
| 480 | + continue; | ||
| 481 | + } | ||
| 482 | + | ||
| 483 | + // Fill event in skymap | ||
| 484 | + map(pixel) += 1.0; | ||
| 485 | + num_in_map++; | ||
| 486 | + | ||
| 487 | + } // endfor: looped over all events | ||
| 488 | + | ||
| 489 | + // Log binning results | ||
| 490 | + if (logTerse()) { | ||
| 491 | + log << std::endl; | ||
| 492 | + log.header1("Mapping"); | ||
| 493 | + log << parformat("Events in list"); | ||
| 494 | + log << obs->events()->size() << std::endl; | ||
| 495 | + log << parformat("Events in map"); | ||
| 496 | + log << num_in_map << std::endl; | ||
| 497 | + log << parformat("Events outside map area"); | ||
| 498 | + log << num_outside_map << std::endl; | ||
| 499 | + log << parformat("Events outside energy range"); | ||
| 500 | + log << num_outside_erange << std::endl; | ||
| 501 | + } | ||
| 502 | + | ||
| 503 | + // Log map | ||
| 504 | + if (logTerse()) { | ||
| 505 | + log << std::endl; | ||
| 506 | + log.header1("Sky map"); | ||
| 507 | + log << map << std::endl; | ||
| 508 | + } | ||
| 509 | + | ||
| 510 | + } // endif: observation was valid | ||
| 511 | + | ||
| 512 | + // Return | ||
| 513 | + return; | ||
| 514 | +} | ||
| 515 | + | ||
| 516 | + | ||
| 517 | +/*========================================================================== | ||
| 518 | + = = | ||
| 519 | + = Private methods = | ||
| 520 | + = = | ||
| 521 | + ==========================================================================*/ | ||
| 522 | + | ||
| 523 | +/***********************************************************************//** | ||
| 524 | + * @brief Initialise class members | ||
| 525 | + ***************************************************************************/ | ||
| 526 | +void ctskymap::init_members(void) | ||
| 527 | +{ | ||
| 528 | + // Initialise members | ||
| 529 | + m_evfile.clear(); | ||
| 530 | + m_outfile.clear(); | ||
| 531 | + m_proj.clear(); | ||
| 532 | + m_coordsys.clear(); | ||
| 533 | + m_obs.clear(); | ||
| 534 | + m_emin = 0.0; | ||
| 535 | + m_emax = 0.0; | ||
| 536 | + m_xref = 9999.0; // Flags unset (use pointing direction) | ||
| 537 | + m_yref = 9999.0; // Flags unset (use pointing direction) | ||
| 538 | + m_binsz = 0.0; | ||
| 539 | + m_nxpix = 0; | ||
| 540 | + m_nypix = 0; | ||
| 541 | + | ||
| 542 | + // Set logger properties | ||
| 543 | + log.date(true); | ||
| 544 | + | ||
| 545 | + // Return | ||
| 546 | + return; | ||
| 547 | +} | ||
| 548 | + | ||
| 549 | + | ||
| 550 | +/***********************************************************************//** | ||
| 551 | + * @brief Copy class members | ||
| 552 | + * | ||
| 553 | + * @param[in] app Application. | ||
| 554 | + ***************************************************************************/ | ||
| 555 | +void ctskymap::copy_members(const ctskymap& app) | ||
| 556 | +{ | ||
| 557 | + // Copy attributes | ||
| 558 | + m_evfile = app.m_evfile; | ||
| 559 | + m_outfile = app.m_outfile; | ||
| 560 | + m_proj = app.m_proj; | ||
| 561 | + m_coordsys = app.m_coordsys; | ||
| 562 | + m_obs = app.m_obs; | ||
| 563 | + m_emin = app.m_emin; | ||
| 564 | + m_emax = app.m_emax; | ||
| 565 | + m_xref = app.m_xref; | ||
| 566 | + m_yref = app.m_yref; | ||
| 567 | + m_binsz = app.m_binsz; | ||
| 568 | + m_nxpix = app.m_nxpix; | ||
| 569 | + m_nypix = app.m_nypix; | ||
| 570 | + | ||
| 571 | + // Return | ||
| 572 | + return; | ||
| 573 | +} | ||
| 574 | + | ||
| 575 | + | ||
| 576 | +/***********************************************************************//** | ||
| 577 | + * @brief Delete class members | ||
| 578 | + ***************************************************************************/ | ||
| 579 | +void ctskymap::free_members(void) | ||
| 580 | +{ | ||
| 581 | + // Write separator into logger | ||
| 582 | + if (logTerse()) | ||
| 583 | + log << std::endl; | ||
| 584 | + | ||
| 585 | + // Return | ||
| 586 | + return; | ||
| 587 | +} |
src/ctskymap/ctskymap.hpp
0 → 100644
| 1 | +/*************************************************************************** | ||
| 2 | + * ctskymap - CTA sky mapping tool * | ||
| 3 | + * ----------------------------------------------------------------------- * | ||
| 4 | + * copyright (C) 2010-2011 by Jurgen Knodlseder * | ||
| 5 | + * ----------------------------------------------------------------------- * | ||
| 6 | + * * | ||
| 7 | + * This program is free software: you can redistribute it and/or modify * | ||
| 8 | + * it under the terms of the GNU General Public License as published by * | ||
| 9 | + * the Free Software Foundation, either version 3 of the License, or * | ||
| 10 | + * (at your option) any later version. * | ||
| 11 | + * * | ||
| 12 | + * This program is distributed in the hope that it will be useful, * | ||
| 13 | + * but WITHOUT ANY WARRANTY; without even the implied warranty of * | ||
| 14 | + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * | ||
| 15 | + * GNU General Public License for more details. * | ||
| 16 | + * * | ||
| 17 | + * You should have received a copy of the GNU General Public License * | ||
| 18 | + * along with this program. If not, see <http://www.gnu.org/licenses/>. * | ||
| 19 | + * * | ||
| 20 | + ***************************************************************************/ | ||
| 21 | +/** | ||
| 22 | + * @file ctskymap.hpp | ||
| 23 | + * @brief CTA sky mapping tool definition | ||
| 24 | + * @author J. Knodlseder | ||
| 25 | + */ | ||
| 26 | + | ||
| 27 | +#ifndef CTSKYMAP_HPP | ||
| 28 | +#define CTSKYMAP_HPP | ||
| 29 | + | ||
| 30 | +/* __ Includes ___________________________________________________________ */ | ||
| 31 | +#include "GammaLib.hpp" | ||
| 32 | +#include "GCTALib.hpp" | ||
| 33 | + | ||
| 34 | +/* __Definitions _________________________________________________________ */ | ||
| 35 | +#define CTSKYMAP_NAME "ctskymap" | ||
| 36 | +#define CTSKYMAP_VERSION "00-01-00" | ||
| 37 | + | ||
| 38 | + | ||
| 39 | +/***********************************************************************//** | ||
| 40 | + * @class ctskymap | ||
| 41 | + * | ||
| 42 | + * @brief CTA sky mapping tool interface defintion | ||
| 43 | + * | ||
| 44 | + * This class creates a sky map from a CTA event list. | ||
| 45 | + ***************************************************************************/ | ||
| 46 | +class ctskymap : public GApplication { | ||
| 47 | +public: | ||
| 48 | + // Constructors and destructors | ||
| 49 | + ctskymap(void); | ||
| 50 | + explicit ctskymap(GObservations obs); | ||
| 51 | + ctskymap(int argc, char *argv[]); | ||
| 52 | + ctskymap(const ctskymap& app); | ||
| 53 | + virtual ~ctskymap(void); | ||
| 54 | + | ||
| 55 | + // Operators | ||
| 56 | + ctskymap& operator= (const ctskymap& app); | ||
| 57 | + | ||
| 58 | + // Methods | ||
| 59 | + void clear(void); | ||
| 60 | + void execute(void); | ||
| 61 | + void run(void); | ||
| 62 | + void save(void); | ||
| 63 | + GObservations& obs(void) { return m_obs; } | ||
| 64 | + void get_parameters(void); | ||
| 65 | + void map_events(GCTAObservation* obs); | ||
| 66 | + | ||
| 67 | +protected: | ||
| 68 | + // Protected methods | ||
| 69 | + void init_members(void); | ||
| 70 | + void copy_members(const ctskymap& app); | ||
| 71 | + void free_members(void); | ||
| 72 | + | ||
| 73 | + // User parameters | ||
| 74 | + std::string m_evfile; //!< Input event list | ||
| 75 | + std::string m_outfile; //!< Output counts map | ||
| 76 | + double m_emin; //!< Minimum energy (TeV) | ||
| 77 | + double m_emax; //!< Maximum energy (TeV) | ||
| 78 | + std::string m_proj; //!< WCS projection | ||
| 79 | + std::string m_coordsys; //!< Coordinate system | ||
| 80 | + double m_xref; //!< Longitude reference coordinate | ||
| 81 | + double m_yref; //!< Latitude reference coordinate | ||
| 82 | + double m_binsz; //!< Pixel size | ||
| 83 | + int m_nxpix; //!< Number of pixels in longitude | ||
| 84 | + int m_nypix; //!< Number of pixels in latitude | ||
| 85 | + | ||
| 86 | + // Protected members | ||
| 87 | + GObservations m_obs; //!< Observation container | ||
| 88 | +}; | ||
| 89 | + | ||
| 90 | +#endif /* CTSKYMAP_HPP */ |
src/ctskymap/ctskymap.par
0 → 100644
| 1 | +############################################################################# | ||
| 2 | +# # | ||
| 3 | +# CTA Science Analysis Tools # | ||
| 4 | +# # | ||
| 5 | +# ------------------------------------------------------------------------- # | ||
| 6 | +# # | ||
| 7 | +# File: ctskymap.par # | ||
| 8 | +# Version: 00-01-00 # | ||
| 9 | +# # | ||
| 10 | +# Author: Jurgen Knodlseder # | ||
| 11 | +# knodlseder@cesr.fr # | ||
| 12 | +# IRAP # | ||
| 13 | +# # | ||
| 14 | +# Purpose: Parameter file for the ctskymap application # | ||
| 15 | +# # | ||
| 16 | +############################################################################# | ||
| 17 | +# | ||
| 18 | +# General parameters | ||
| 19 | +#=================== | ||
| 20 | +evfile, f, a, "test/data/crab_events.fits.gz",,,"Event data file name" | ||
| 21 | +outfile, f, a, "skymap.fits",,, "Output file name" | ||
| 22 | +# | ||
| 23 | +# Event selection parameters | ||
| 24 | +#=========================== | ||
| 25 | +emin, r, a, 0.1,,,"Minimum energy in TeV" | ||
| 26 | +emax, r, a, 100.0,,,"Maximum energy in TeV" | ||
| 27 | +# | ||
| 28 | +# Spatial binning parameters | ||
| 29 | +#=========================== | ||
| 30 | +nxpix, i, a, 200,,, "Size of the X axis in pixels" | ||
| 31 | +nypix, i, a, 200,,, "Size of the Y axis in pixels" | ||
| 32 | +binsz, r, a, 0.02,,, "Image scale (in degrees/pixel)" | ||
| 33 | +coordsys, s, a, "CEL", CEL|GAL,,"Coordinate system (CEL - celestial, GAL - galactic)" | ||
| 34 | +xref, r, a, 83.63,,, "First coordinate of image center in degrees (RA or galactic l)" | ||
| 35 | +yref, r, a, 22.01,,, "Second coordinate of image center in degrees (DEC or galactic b)" | ||
| 36 | +axisrot, r, h, 0.0,,, "Rotation angle of image axis, in degrees" | ||
| 37 | +proj, s, a, "CAR",,, "Projection method e.g. AIT|ARC|CAR|GLS|MER|NCP|SIN|STG|TAN" | ||
| 38 | + | ||
| 39 | +# | ||
| 40 | +# Standard parameters | ||
| 41 | +#==================== | ||
| 42 | +chatter, i, h, 1,0,4, "Chattiness of output" | ||
| 43 | +clobber, b, h, yes,,, "Overwrite existing output files with new output files?" | ||
| 44 | +debug, b, h, no,,, "Debugging mode activated" | ||
| 45 | +mode, s, h, "ql",,, "Mode of automatic parameters" |
src/ctskymap/main.cpp
0 → 100644
| 1 | +/*************************************************************************** | ||
| 2 | + * ctskymap - CTA sky mapping main code * | ||
| 3 | + * ----------------------------------------------------------------------- * | ||
| 4 | + * copyright (C) 2010-2011 by Jurgen Knodlseder * | ||
| 5 | + * ----------------------------------------------------------------------- * | ||
| 6 | + * * | ||
| 7 | + * This program is free software: you can redistribute it and/or modify * | ||
| 8 | + * it under the terms of the GNU General Public License as published by * | ||
| 9 | + * the Free Software Foundation, either version 3 of the License, or * | ||
| 10 | + * (at your option) any later version. * | ||
| 11 | + * * | ||
| 12 | + * This program is distributed in the hope that it will be useful, * | ||
| 13 | + * but WITHOUT ANY WARRANTY; without even the implied warranty of * | ||
| 14 | + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * | ||
| 15 | + * GNU General Public License for more details. * | ||
| 16 | + * * | ||
| 17 | + * You should have received a copy of the GNU General Public License * | ||
| 18 | + * along with this program. If not, see <http://www.gnu.org/licenses/>. * | ||
| 19 | + * * | ||
| 20 | + ***************************************************************************/ | ||
| 21 | +/** | ||
| 22 | + * @file main.cpp | ||
| 23 | + * @brief CTA sky mapping tool main code | ||
| 24 | + * @author J. Knodlseder | ||
| 25 | + */ | ||
| 26 | + | ||
| 27 | +/* __ Includes ___________________________________________________________ */ | ||
| 28 | +#ifdef HAVE_CONFIG_H | ||
| 29 | +#include <config.h> | ||
| 30 | +#endif | ||
| 31 | +#include <cstdio> | ||
| 32 | +#include "ctskymap.hpp" | ||
| 33 | + | ||
| 34 | + | ||
| 35 | +/***********************************************************************//** | ||
| 36 | + * @brief Main entry point of application | ||
| 37 | + * | ||
| 38 | + * @param[in] argc Number of command line arguments. | ||
| 39 | + * @param[in] argv Command line arguments. | ||
| 40 | + * | ||
| 41 | + * This is the main entry point of the ctskymap application. It allocates a | ||
| 42 | + * ctskymap object and runs the application. | ||
| 43 | + ***************************************************************************/ | ||
| 44 | +int main (int argc, char *argv[]) | ||
| 45 | +{ | ||
| 46 | + // Initialise return code | ||
| 47 | + int rc = 1; | ||
| 48 | + | ||
| 49 | + // Create instance of application | ||
| 50 | + ctskymap application(argc, argv); | ||
| 51 | + | ||
| 52 | + // Run application | ||
| 53 | + try { | ||
| 54 | + // Execute application | ||
| 55 | + application.execute(); | ||
| 56 | + | ||
| 57 | + // Signal success | ||
| 58 | + rc = 0; | ||
| 59 | + } | ||
| 60 | + catch (std::exception &e) { | ||
| 61 | + | ||
| 62 | + // Extract error message | ||
| 63 | + std::string message = e.what(); | ||
| 64 | + std::string signal = "*** ERROR encounterted in the execution of" | ||
| 65 | + " ctskymap. Run aborted ..."; | ||
| 66 | + | ||
| 67 | + // Write error in logger | ||
| 68 | + application.log << signal << std::endl; | ||
| 69 | + application.log << message << std::endl; | ||
| 70 | + | ||
| 71 | + // Write error on standard output | ||
| 72 | + std::cout << signal << std::endl; | ||
| 73 | + std::cout << message << std::endl; | ||
| 74 | + | ||
| 75 | + } // endcatch: catched any application error | ||
| 76 | + | ||
| 77 | + // Return | ||
| 78 | + return rc; | ||
| 79 | +} |