Commit f349fd81b4ec7d9625aeaaa14f8bd945cfebc838

Authored by Simon Bonnefoy
Committed by Jürgen Knödlseder
1 parent 93f8b5d8

fix the get_variability_sig in ctfindvar and the src/Makefile.am (#2713)

src/Makefile.am
@@ -60,7 +60,6 @@ pkginclude_HEADERS = ctools.hpp @@ -60,7 +60,6 @@ pkginclude_HEADERS = ctools.hpp
60 lib_LTLIBRARIES = libctools.la 60 lib_LTLIBRARIES = libctools.la
61 libctools_la_SOURCES = 61 libctools_la_SOURCES =
62 libctools_la_LIBADD = support/libsupport.la \ 62 libctools_la_LIBADD = support/libsupport.la \
63 - ctfindvar/libctfindvar.la \  
64 ctbin/libctbin.la \ 63 ctbin/libctbin.la \
65 ctlike/libctlike.la \ 64 ctlike/libctlike.la \
66 ctmodel/libctmodel.la \ 65 ctmodel/libctmodel.la \
src/ctfindvar/ctfindvar.cpp
@@ -31,6 +31,7 @@ @@ -31,6 +31,7 @@
31 #include "ctfindvar.hpp" 31 #include "ctfindvar.hpp"
32 #include "GammaLib.hpp" 32 #include "GammaLib.hpp"
33 #include "GCTALib.hpp" 33 #include "GCTALib.hpp"
  34 +#include <cmath>
34 35
35 /* __ OpenMP section _____________________________________________________ */ 36 /* __ OpenMP section _____________________________________________________ */
36 //#ifdef _OPENMP 37 //#ifdef _OPENMP
@@ -233,7 +234,7 @@ void ctfindvar::run(void) @@ -233,7 +234,7 @@ void ctfindvar::run(void)
233 // GNdarray GNdarray_temp(NBINS); 234 // GNdarray GNdarray_temp(NBINS);
234 // for (int i=0, i<number of pixels; i++) 235 // for (int i=0, i<number of pixels; i++)
235 // { 236 // {
236 - // get_variability_sig(pix_number, GNdarray_temp) 237 + // get_variability_sig(pix_number,nbins, G Ndarray_temp)
237 // if (position == source position) 238 // if (position == source position)
238 // { 239 // {
239 // GNdarray_source= GNdarray_temp; 240 // GNdarray_source= GNdarray_temp;
@@ -254,7 +255,9 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number, const int&amp; nbins, GNd @@ -254,7 +255,9 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number, const int&amp; nbins, GNd
254 { 255 {
255 //const int nbins = hist->GetNbinsX(); 256 //const int nbins = hist->GetNbinsX();
256 //const int nbins = VECTOR.size() 257 //const int nbins = VECTOR.size()
257 - //const int nbins = GSKYMAP.maps(); 258 + //const int nbins = m_counts.maps();
  259 + std::vector<bool> accepted_bin_bckg_vector;
  260 + std::vector<double> excess_bin_vector;
258 int background_bin_array[nbins]; 261 int background_bin_array[nbins];
259 bool background_validated=false; 262 bool background_validated=false;
260 int non, noff; 263 int non, noff;
@@ -283,7 +286,7 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number, const int&amp; nbins, GNd @@ -283,7 +286,7 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number, const int&amp; nbins, GNd
283 { 286 {
284 //background_count+=hist->GetBinContent(j+1); 287 //background_count+=hist->GetBinContent(j+1);
285 //background_count+=VECTOR(j+1).first(pix_number); 288 //background_count+=VECTOR(j+1).first(pix_number);
286 - background_count+=GSKYMAP(pix_number, j); 289 + background_count+= m_counts(pix_number, j);
287 alpha++; 290 alpha++;
288 } 291 }
289 } 292 }
@@ -291,10 +294,20 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number, const int&amp; nbins, GNd @@ -291,10 +294,20 @@ void ctfindvar::get_variability_sig(const int&amp; pix_number, const int&amp; nbins, GNd
291 background_bin_array[i] = background_count; //The background is averaged on the number of bins -1 294 background_bin_array[i] = background_count; //The background is averaged on the number of bins -1
292 //non = hist->GetBinContent(i+1); 295 //non = hist->GetBinContent(i+1);
293 //non = VECTOR(i+1).first(pix_number); 296 //non = VECTOR(i+1).first(pix_number);
294 - non = GSKYMAP(pix_number, i); 297 + non = m_counts(pix_number, i);
295 noff = background_bin_array[i]; 298 noff = background_bin_array[i];
296 alpha = (1./alpha); 299 alpha = (1./alpha);
297 - sig = Utilities::Statistics::Significance(non,noff,alpha); 300 +
  301 + /////////////////////////////////////////////////////////////////////////////////
  302 + // Compute sensitivity in Gaussian sigma
  303 + double alpha1 = alpha + 1.0;
  304 + double ntotal = non+noff;
  305 + double arg1 = non/ntotal;
  306 + double arg2 = noff/ntotal;
  307 + double term1 = non * std::log((alpha1/alpha)*arg1);
  308 + double term2 = noff * std::log(alpha1*arg2);
  309 + sig = sqrt(2.0 * (term1 + term2));
  310 + /////////////////////////////////////////////////
298 //sig_bin_vector[i]=sig; 311 //sig_bin_vector[i]=sig;
299 sig_histogram(i)=sig; 312 sig_histogram(i)=sig;
300 excess_bin_vector[i]=non - alpha*noff; 313 excess_bin_vector[i]=non - alpha*noff;
@@ -451,4 +464,4 @@ void ctfindvar::init_gtis(void) @@ -451,4 +464,4 @@ void ctfindvar::init_gtis(void)
451 464
452 // ToDO: This is wrong only making 1 bin 465 // ToDO: This is wrong only making 1 bin
453 m_gti = GGti(tstart, tstop); 466 m_gti = GGti(tstart, tstop);
454 -}  
455 \ No newline at end of file 467 \ No newline at end of file
  468 +}