Commit 5aa6adb3cbc49ef1a91cb9932bf6b45fd8500f44

Authored by Jürgen Knödlseder
1 parent a1d911f2

Integrate code

src/ctselect/ctselect.cpp
@@ -478,7 +478,7 @@ void ctselect::init_members(void) @@ -478,7 +478,7 @@ void ctselect::init_members(void)
478 m_phases.clear(); 478 m_phases.clear();
479 m_select_energy = false; 479 m_select_energy = false;
480 m_select_phase = false; 480 m_select_phase = false;
481 - m_forcesel = false; 481 + m_forcesel = false;
482 482
483 // Return 483 // Return
484 return; 484 return;
@@ -507,7 +507,7 @@ void ctselect::copy_members(const ctselect& app) @@ -507,7 +507,7 @@ void ctselect::copy_members(const ctselect& app)
507 m_phases = app.m_phases; 507 m_phases = app.m_phases;
508 m_select_energy = app.m_select_energy; 508 m_select_energy = app.m_select_energy;
509 m_select_phase = app.m_select_phase; 509 m_select_phase = app.m_select_phase;
510 - m_forcesel = app.m_forcesel; 510 + m_forcesel = app.m_forcesel;
511 511
512 // Return 512 // Return
513 return; 513 return;
@@ -550,9 +550,10 @@ void ctselect::get_parameters(void) @@ -550,9 +550,10 @@ void ctselect::get_parameters(void)
550 // Query the RoI 550 // Query the RoI
551 GCTAPointing pnt; 551 GCTAPointing pnt;
552 GCTARoi roi = get_roi(pnt); 552 GCTARoi roi = get_roi(pnt);
553 - // Signal whether the RoI selection should be enforced  
554 - // even when leading to invalid RoI  
555 - m_forcesel = (has_par("forcesel") && (*this)["forcesel"].boolean()); 553 +
  554 + // Signal whether the RoI selection should be enforced even when leading
  555 + // to invalid RoI
  556 + m_forcesel = (*this)["forcesel"].boolean();
556 557
557 // Query Good Time Intervals. The GammaLib time reference is specified 558 // Query Good Time Intervals. The GammaLib time reference is specified
558 // as a dummy argument since the relevant intervals will be queried later 559 // as a dummy argument since the relevant intervals will be queried later
@@ -735,7 +736,8 @@ void ctselect::select_events(GCTAObservation* obs, @@ -735,7 +736,8 @@ void ctselect::select_events(GCTAObservation* obs,
735 // Signal that RoI selection should be performed 736 // Signal that RoI selection should be performed
736 bool select_roi = roi.is_valid(); 737 bool select_roi = roi.is_valid();
737 738
738 - //Create boolean to remember if RoI selection is enforced regardless of its validity 739 + // Create boolean to remember if RoI selection is enforced regardless of
  740 + // its validity
739 bool enforced_roi = false; 741 bool enforced_roi = false;
740 742
741 // Set time selection interval. We make sure here that the time selection 743 // Set time selection interval. We make sure here that the time selection
@@ -921,37 +923,44 @@ void ctselect::select_events(GCTAObservation* obs, @@ -921,37 +923,44 @@ void ctselect::select_events(GCTAObservation* obs,
921 GSkyDir centre; 923 GSkyDir centre;
922 centre.radec_deg(ra, dec); 924 centre.radec_deg(ra, dec);
923 double distance = centre.dist_deg(roi_centre); 925 double distance = centre.dist_deg(roi_centre);
924 - // If the requested RoI and selected RoI are centred at the same position  
925 - // and the requested RoI is bigger then the original 926 +
  927 + // If the requested RoI and selected RoI are centred at the same
  928 + // position and the requested RoI is bigger then the original
926 // keep original radius 929 // keep original radius
927 // RoIs same centre = distance < 1.e-4 930 // RoIs same centre = distance < 1.e-4
928 // to cope with numerical precision limitations 931 // to cope with numerical precision limitations
929 - if (distance < 1.e-4 && rad > roi_radius) { 932 + if ((distance < 1.e-4) && (rad > roi_radius)) {
930 rad = roi_radius; 933 rad = roi_radius;
931 } 934 }
932 - //Otherwise, check if requested RoI is enclosed within original RoI  
933 - else {  
934 - if (distance + rad > roi_radius) {  
935 - //If selection is enforced log warning and store  
936 - if (m_forcesel){ 935 +
  936 + // ... otherwise, check if requested RoI is enclosed within
  937 + // original RoI
  938 + else if ((distance + rad) > roi_radius) {
  939 +
  940 + // If selection is enforced log warning and store
  941 + if (m_forcesel) {
937 enforced_roi = true; 942 enforced_roi = true;
938 - log_string(TERSE, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");  
939 - log_string(TERSE,"Enforced RoI selection leading to inconsistent RoI. "\  
940 - "RoI of data: Centre(RA,DEC)=("+gammalib::str(ra)+", "+  
941 - gammalib::str(dec)+") deg, Radius="+gammalib::str(rad)+  
942 - " deg. Requested RoI: Centre(RA,DEC)=("+  
943 - gammalib::str(roi_centre.ra_deg())+", "+  
944 - gammalib::str(roi_centre.dec_deg())+") deg, Radius="+  
945 - gammalib::str(roi_radius)+" deg. The resulting event list "\  
946 - "file is not usable for likelihood analysis");  
947 - log_string(TERSE, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");  
948 - }  
949 - // Otherwise, throw an exception  
950 - else {  
951 - std::string message;  
952 - message = "Invalid RoI selection: the new RoI must be enclosed in the original RoI";  
953 - throw GException::invalid_value(G_RUN, message);  
954 - } 943 + std::string msg = \
  944 + "WARNING: Enforced RoI selection leading to inconsis" \
  945 + "tent RoI.\n" \
  946 + " RoI of data: Centre(RA,DEC)=("+ \
  947 + gammalib::str(ra)+", "+gammalib::str(dec)+") deg, " \
  948 + "Radius="+gammalib::str(rad)+" deg.\n" \
  949 + " Requested RoI: Centre(RA,DEC)=("+ \
  950 + gammalib::str(roi_centre.ra_deg())+", "+ \
  951 + gammalib::str(roi_centre.dec_deg())+") deg, " \
  952 + "Radius="+gammalib::str(roi_radius)+" deg.\n" \
  953 + " The resulting event list file is not usable " \
  954 + "for likelihood analysis";
  955 + log_string(TERSE, msg);
  956 + }
  957 +
  958 + // ... otherwise, throw an exception
  959 + else {
  960 + std::string msg = \
  961 + "Invalid RoI selection: the new RoI must be enclosed "
  962 + "in the original RoI";
  963 + throw GException::invalid_value(G_RUN, msg);
955 } 964 }
956 } 965 }
957 } 966 }
@@ -971,7 +980,7 @@ void ctselect::select_events(GCTAObservation* obs, @@ -971,7 +980,7 @@ void ctselect::select_events(GCTAObservation* obs,
971 gammalib::str(ra)+", "+gammalib::str(dec)+") deg, " 980 gammalib::str(ra)+", "+gammalib::str(dec)+") deg, "
972 "Radius="+gammalib::str(rad)+" deg"); 981 "Radius="+gammalib::str(rad)+" deg");
973 } 982 }
974 - 983 +
975 // Format RoI selection 984 // Format RoI selection
976 char cra[80]; 985 char cra[80];
977 char cdec[80]; 986 char cdec[80];
@@ -1049,7 +1058,8 @@ void ctselect::select_events(GCTAObservation* obs, @@ -1049,7 +1058,8 @@ void ctselect::select_events(GCTAObservation* obs,
1049 // closed later 1058 // closed later
1050 list->fetch(); 1059 list->fetch();
1051 1060
1052 - // If RoI selection has been applied without forcing then set the event list RoI 1061 + // If RoI selection has been applied without forcing then set the event
  1062 + // list RoI
1053 if (select_roi && !enforced_roi) { 1063 if (select_roi && !enforced_roi) {
1054 GCTAInstDir instdir; 1064 GCTAInstDir instdir;
1055 instdir.dir().radec_deg(ra, dec); 1065 instdir.dir().radec_deg(ra, dec);
src/ctselect/ctselect.hpp
@@ -81,7 +81,6 @@ protected: @@ -81,7 +81,6 @@ protected:
81 81
82 // User parameters 82 // User parameters
83 std::string m_outobs; //!< Output event list or XML file 83 std::string m_outobs; //!< Output event list or XML file
84 - //GGti m_gti; //!< Time selection  
85 double m_emin; //!< Lower energy 84 double m_emin; //!< Lower energy
86 double m_emax; //!< Upper energy 85 double m_emax; //!< Upper energy
87 std::string m_expr; //!< Selection expression 86 std::string m_expr; //!< Selection expression
@@ -93,8 +92,6 @@ protected: @@ -93,8 +92,6 @@ protected:
93 std::vector<std::string> m_infiles; //!< Input event filenames 92 std::vector<std::string> m_infiles; //!< Input event filenames
94 std::vector<std::string> m_evtname; //!< Event extension names 93 std::vector<std::string> m_evtname; //!< Event extension names
95 std::vector<std::string> m_gtiname; //!< GTI extension names 94 std::vector<std::string> m_gtiname; //!< GTI extension names
96 - //GTime m_timemin; //!< Earliest time  
97 - //GTime m_timemax; //!< Latest time  
98 GPhases m_phases; //!< Phase intervals 95 GPhases m_phases; //!< Phase intervals
99 bool m_select_energy; //!< Perform energy selection 96 bool m_select_energy; //!< Perform energy selection
100 bool m_select_phase; //!< Perform phase selection 97 bool m_select_phase; //!< Perform phase selection
src/ctselect/ctselect.par
@@ -23,18 +23,18 @@ prefix, s, h, selected_,,, &quot;Prefix for event lists in observation definition XML @@ -23,18 +23,18 @@ prefix, s, h, selected_,,, &quot;Prefix for event lists in observation definition XML
23 # 23 #
24 # Selection parameters 24 # Selection parameters
25 #===================== 25 #=====================
26 -usepnt, b, h, yes,,, "Use pointing instead of RA/DEC parameters?"  
27 -ra, r, a, 83.63,0,360, "RA for ROI centre (degrees)"  
28 -dec, r, a, 22.51,-90,90, "Dec for ROI centre (degrees)"  
29 -rad, r, a, 3.0,0,180, "Radius of ROI around pointing or specified RA/DEC (degrees)"  
30 -forcesel,b, h, no,,, "Force ROI selection even if it produces invalid ROI?"  
31 -tmin, t, a, NONE,,, "Start time (UTC string, JD, MJD or MET in seconds)"  
32 -tmax, t, a, NONE,,, "Stop time (UTC string, JD, MJD or MET in seconds)"  
33 -emin, r, a, 0.1,,, "Lower energy limit (TeV)"  
34 -emax, r, a, 100.0,,, "Upper energy limit (TeV)"  
35 -phase, s, h, NONE,,, "Phase expression in the format phasemin0:phasemax0,phasemin1:phasemax1,..."  
36 -expr, s, h, ,,, "Additional selection expression"  
37 -usethres,s, h, NONE,NONE|DEFAULT|USER,, "Energy threshold type" 26 +usepnt, b, h, yes,,, "Use pointing instead of RA/DEC parameters?"
  27 +ra, r, a, 83.63,0,360, "RA for ROI centre (degrees)"
  28 +dec, r, a, 22.51,-90,90, "Dec for ROI centre (degrees)"
  29 +rad, r, a, 3.0,0,180, "Radius of ROI around pointing or specified RA/DEC (degrees)"
  30 +forcesel, b, h, no,,, "Force ROI selection even if it produces invalid ROI?"
  31 +tmin, t, a, NONE,,, "Start time (UTC string, JD, MJD or MET in seconds)"
  32 +tmax, t, a, NONE,,, "Stop time (UTC string, JD, MJD or MET in seconds)"
  33 +emin, r, a, 0.1,,, "Lower energy limit (TeV)"
  34 +emax, r, a, 100.0,,, "Upper energy limit (TeV)"
  35 +phase, s, h, NONE,,, "Phase expression in the format phasemin0:phasemax0,phasemin1:phasemax1,..."
  36 +expr, s, h, ,,, "Additional selection expression"
  37 +usethres, s, h, NONE,NONE|DEFAULT|USER,, "Energy threshold type"
38 38
39 # 39 #
40 # Standard parameters 40 # Standard parameters