Karim Sabri / ctools

Commit 056df92382d28a96e45df588cb02227b50e5c911

Authored by Josh Cardenzana
Committed by Jürgen Knödlseder
1 parent 29669ee5

Apply style formatting according to ctools standard (#2713)

Inline Side-by-side
Showing 2 changed files with 61 additions and 66 deletions
src/ctfindvar/ctfindvar.cpp
  View file @056df92
... ... @@ -32,7 +32,6 @@
32 32 #include "GammaLib.hpp"
33 33 #include "GCTALib.hpp"
34 34 #include <cmath>
35   -#include <sstream>
36 35  
37 36 /* __ OpenMP section _____________________________________________________ */
38 37 #ifdef _OPENMP
... ... @@ -288,12 +287,9 @@ void ctfindvar::run(void)
288 287 // ... otherwise only one source position is queried
289 288 else {
290 289 GSkyDir srcSkyDir;
291   - if ((*this)["coordsys"].string()=="CEL")
292   - {
  290 + if ((*this)["coordsys"].string() == "CEL") {
293 291 srcSkyDir.radec_deg((*this)["xsrc"].real(),(*this)["ysrc"].real());
294   - }
295   - else
296   - {
  292 + } else {
297 293 srcSkyDir.lb_deg((*this)["xsrc"].real(),(*this)["ysrc"].real());
298 294 }
299 295  
... ... @@ -308,13 +304,13 @@ void ctfindvar::run(void)
308 304 m_pixsigmax = GNdarray(nbins);
309 305  
310 306 //Prepare the final skymap with the max significance of each pixel
311   - double max_sig=0;
  307 + double max_sig = 0;
312 308  
313 309 log_header1(NORMAL, "Looping Over Pixels");
314 310  
315 311 //looping over all the pixels in the cube
316 312 #pragma omp parallel for
317   - for (int pix_number=0; pix_number<m_counts.npix(); pix_number++) {
  313 + for (int pix_number = 0; pix_number < m_counts.npix(); pix_number++) {
318 314  
319 315 //Preparing array to store significance evolution of
320 316 //each pixel
... ... @@ -322,7 +318,7 @@ void ctfindvar::run(void)
322 318 double total_counts=0;
323 319  
324 320 //Checking the integrated count of each pixel
325   - for (int k=0;k<m_counts.nmaps();k++) {
  321 + for (int k = 0; k < m_counts.nmaps(); k++) {
326 322 total_counts += m_counts(pix_number, k);
327 323 }
328 324 #ifdef G_DEBUG
... ... @@ -337,13 +333,13 @@ void ctfindvar::run(void)
337 333 get_variability_sig(pix_number,nbins, pixSig);
338 334  
339 335 //Getting the significance evolution for the source
340   - for (int src=0; src<srcInxPix.size(); src++) {
  336 + for (int src = 0; src < srcInxPix.size(); src++) {
341 337  
342 338 // Store the distribution if the source is located at this position
343 339 if (srcInxPix[src] == pix_number) {
344 340  
345 341 #pragma omp critical(ctfind_run)
346   - for (int i=0; i<nbins; i++) {
  342 + for (int i = 0; i < nbins; i++) {
347 343 m_pixsigsrc(src,i) = pixSig(i);
348 344 }
349 345 #ifdef G_DEBUG
... ... @@ -371,8 +367,7 @@ void ctfindvar::run(void)
371 367 // Getting the evolution for the pix with highest significance
372 368 #pragma omp critical(ctfind_run)
373 369 {
374   - if (max(pixSig) > max_sig)
375   - {
  370 + if (max(pixSig) > max_sig) {
376 371 // Update the pixel position with the maximum sigma
377 372 max_sig = max(pixSig);
378 373 m_max_sig_dir = m_counts.inx2dir(pix_number);
... ... @@ -384,7 +379,7 @@ void ctfindvar::run(void)
384 379  
385 380 //If the significance is greater than sig threshold, the position
386 381 //is stored in a model
387   - if (max(pixSig)> m_sig_threshold) {
  382 + if (max(pixSig) > m_sig_threshold) {
388 383 GSkyDir dir_pix = m_counts.inx2dir(pix_number);
389 384 fill_model_sig_pos(dir_pix);
390 385 }
... ... @@ -411,9 +406,8 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number,
411 406 const int& nbins,
412 407 GNdarray& sig_histogram)
413 408 {
414   - std::vector<bool> accepted_bin_bckg_vector(nbins, true);
415   - std::vector<double> excess_bin_vector(nbins, 0.0);
416   - bool background_validated=false;
  409 + std::vector<bool> accepted_bin_bckg_vector(nbins, true);
  410 + bool background_validated = false;
417 411 double non, noff;
418 412 double alpha, sig;
419 413  
... ... @@ -421,39 +415,38 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number,
421 415 std::vector<double> alpha_vector(nbins, 0.0);
422 416 fill_alpha_vector(pix_number, alpha_vector);
423 417  
424   - while (background_validated==false)
425   - {
426   - background_validated=true;
  418 + // Loop over pixels until all background pixls are removed
  419 + while (background_validated == false) {
  420 + background_validated = true;
  421 +
  422 + // Loop over all the GTIs of the pixel
  423 + for (int i = 0; i < nbins; i++) {
427 424  
428   - for (int i=0; i< nbins; i++) //looping over all the GTIs of the pixel
429   - {
430   - //the GTI is discared from bckg calculation and not checked again.
431   - if (accepted_bin_bckg_vector[i]==0) {
  425 + // The GTI is discared from bckg calculation and not checked again.
  426 + if (!accepted_bin_bckg_vector[i]) {
432 427 continue;
433 428 }
434 429 // Check if bin fails minoff check or no observations overlap it
435   - else if (m_counts(pix_number, i) < m_minoff || alpha_vector[i]==0) {
436   - accepted_bin_bckg_vector[i]=0;
  430 + else if (m_counts(pix_number, i) < m_minoff || alpha_vector[i] == 0) {
  431 + accepted_bin_bckg_vector[i] = false;
437 432 continue;
438 433 }
439 434  
440 435 noff = 0.0;
441   - alpha = 0;
442   - for (int j=0;j<nbins;j++) // for one GTI selected (i), looping over all the others (j).
443   - {
444   - if (j!=i && accepted_bin_bckg_vector[j]==1)
445   - {
  436 + alpha = 0.0;
  437 + // Gor one GTI selected (i), looping over all the others (j).
  438 + for (int j = 0; j < nbins; j++) {
  439 + if (j != i && accepted_bin_bckg_vector[j] == 1) {
446 440 noff += m_counts(pix_number, j);
447 441 alpha += alpha_vector[j];
448 442 }
449 443 }
450 444  
451   - //The background is averaged on the number of bins -1
  445 + // The background is averaged on the number of bins -1
452 446 non = m_counts(pix_number, i);
453 447 alpha = alpha_vector[i]/alpha;
454 448  
455   - /////////////////////////////////////////////////////////////////////////////////
456   - // Compute sensitivity in Gaussian sigma
  449 + // Compute sensitivity in Gaussian sigma using Li & Ma
457 450 double alpha1 = alpha + 1.0;
458 451 double ntotal = non+noff;
459 452 double arg1 = non/ntotal;
... ... @@ -474,17 +467,16 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number,
474 467 if (alpha == 0.0) {
475 468 sig = 0.0;
476 469 }
477   - /////////////////////////////////////////////////
478   - //sig_bin_vector[i]=sig;
479   - sig_histogram(i)=sig;
480   - excess_bin_vector[i]=non - alpha*noff;
  470 +
  471 + // Update the significance
  472 + sig_histogram(i) = sig;
481 473 #ifdef G_DEBUG
482 474 //std::cout << "significance of the bin : " << i << " : " << sig << " - alpha : " << alpha << " - non: " << non << " - noff: " << noff << "- excess: " << non - alpha*noff << std::endl;
483 475 #endif
484   - if (sig>m_sig_threshold ) // if the bin is significant, it is removed from the bckg and we loop again.
485   - {
486   - accepted_bin_bckg_vector[i]=0;
487   - background_validated=false;
  476 + // If the bin is significant, it is removed from the bckg and we loop again.
  477 + if (sig > m_sig_threshold) {
  478 + accepted_bin_bckg_vector[i] = false;
  479 + background_validated = false;
488 480 }
489 481 }
490 482  
... ... @@ -511,9 +503,9 @@ void ctfindvar::fill_model_sig_pos(const GSkyDir&amp; dir_pix)
511 503 //Extracting the coordinate to name the model
512 504 double ra = dir_pix.ra_deg();
513 505 double dec = dir_pix.dec_deg();
514   - std::stringstream name;
515   - name << "Sky_direction_"<< ra << "_" << dec ;
516   - model.name(name.str());
  506 + std::string name = "Sky_direction_" + gammalib::str(ra, 4) +
  507 + "_" + gammalib::str(dec, 4);
  508 + model.name(name);
517 509  
518 510 //Add the model for this pixel to the output mode
519 511 m_model_above_thr.append(model);
... ... @@ -532,7 +524,7 @@ void ctfindvar::fill_alpha_vector(const int&amp; pix_number,
532 524 GSkyDir pix_dir = m_counts.inx2dir(pix_number);
533 525  
534 526 // Loop over all observations
535   - for (int i=0; i<m_obs.size(); i++) {
  527 + for (int i = 0; i < m_obs.size(); i++) {
536 528 // Skip if observation does not overlap with this pixel position
537 529 GCTAObservation* obs = dynamic_cast<GCTAObservation*>(m_obs[i]);
538 530 GCTARoi roi = obs->roi();
... ... @@ -549,7 +541,7 @@ void ctfindvar::fill_alpha_vector(const int&amp; pix_number,
549 541 const GGti& obs_gti(obs->gti());
550 542  
551 543 // Loop over all time bins
552   - for (int j=0; j<m_gti.size(); j++) {
  544 + for (int j = 0; j < m_gti.size(); j++) {
553 545  
554 546 // Make sure observation overlaps with this time interval
555 547 double exposure = gti_overlap(m_gti[j], obs_gti);
... ... @@ -605,7 +597,7 @@ int ctfindvar::time2inx(const GTime&amp; time)
605 597 int map_index = -1;
606 598  
607 599 // Loop over all GTIs
608   - for (int i=0; i<m_gti.size(); i++) {
  600 + for (int i = 0; i < m_gti.size(); i++) {
609 601 // Check if interval contains the time
610 602 if (m_gti[i].contains(time)) {
611 603 map_index = i;
... ... @@ -665,10 +657,11 @@ void ctfindvar::save(void)
665 657  
666 658 // Write the output file
667 659 GFilename outfilename = prefix + "signifmap.fits";
  660 + log_value(TERSE, "Saving results to ", outfilename);
668 661 m_outfile.saveto(outfilename, (*this)["clobber"].boolean());
669 662  
670 663 //Save the model file with position with significance
671   - if (m_model_above_thr.size() > 2){
  664 + if (m_model_above_thr.size() > 2) {
672 665 GFilename model_filename = prefix + "model_significant_directions.xml";
673 666 m_model_above_thr.save(model_filename);
674 667 }
... ... @@ -987,10 +980,9 @@ void ctfindvar::init_gtis(void)
987 980 throw GException::invalid_value(G_INIT_GTIS,
988 981 "Start time is equal to stop time");
989 982 } else if (bins < 2) {
990   - std::stringstream msg;
991   - msg << "Method requires at least two time bins ("
992   - << bins << " bins found).";
993   - throw GException::invalid_value(G_INIT_GTIS, msg.str());
  983 + std::string msg = "Method requires at least two time bins (" +
  984 + gammalib::str(bins) + " bins found).";
  985 + throw GException::invalid_value(G_INIT_GTIS, msg);
994 986 }
995 987  
996 988 for (int i=0; i<bins; i++) {
... ...
src/ctfindvar/ctfindvar.hpp
  View file @056df92
... ... @@ -68,18 +68,21 @@ public:
68 68  
69 69 protected:
70 70 // Protected methods
71   - void init_cube(void);
72   - void fill_model_sig_pos(const GSkyDir& dir_pix);
73   - void init_gtis(void);
74   - void init_members(void);
75   - void copy_members(const ctfindvar& app);
76   - void free_members(void);
77   - void fill_cube(GCTAObservation* obs);
78   - void fill_alpha_vector(const int& pix_number,
79   - std::vector<double>& alpha_vector);
80   - void get_parameters(void);
81   - void get_variability_sig(const int& pix_number, const int& nbins, GNdarray& sig_histogram);
82   - double gti_overlap(const GGti& gti1, const GGti& gti2);
  71 + void init_cube(void);
  72 + void fill_model_sig_pos(const GSkyDir& dir_pix);
  73 + void init_gtis(void);
  74 + void init_members(void);
  75 + void copy_members(const ctfindvar& app);
  76 + void free_members(void);
  77 + void fill_cube(GCTAObservation* obs);
  78 + void fill_alpha_vector(const int& pix_number,
  79 + std::vector<double>& alpha_vector);
  80 + void get_parameters(void);
  81 + void get_variability_sig(const int& pix_number,
  82 + const int& nbins,
  83 + GNdarray& sig_histogram);
  84 + double gti_overlap(const GGti& gti1,
  85 + const GGti& gti2);
83 86  
84 87 // Writing methods for individual source histograms
85 88 void write_srchist(void);
... ... @@ -97,7 +100,7 @@ protected:
97 100 GSkyDir m_max_sig_dir; //!< Sky direction associated with maximum significance
98 101 double m_minoff; //!< Minimum counts for use in significance calculation
99 102 double m_sig_threshold; //!< Minimum significance required to set source as variable
100   - GFits m_outfile; //!< FITS object to store output
  103 + GFits m_outfile; //!< FITS object to store output
101 104 GSkyMap m_peaksigmap; //!< Skymap holding the maximum significance
102 105 GNdarray m_pixsigsrc; //!< Store distributions of the source significances
103 106 GNdarray m_pixsigmax; //!< Store distribution for pixel with max significance
... ...