Commit 7a6df46bd8e0a1493747e1f9747ead2b3c05743b
1 parent
df35383b
add require_inobs for observation setup in csresmap
Showing
2 changed files
with
59 additions
and
49 deletions
scripts/cspull.py
| ... | ... | @@ -19,6 +19,7 @@ |
| 19 | 19 | # |
| 20 | 20 | # ========================================================================== |
| 21 | 21 | import gammalib |
| 22 | +import ctools | |
| 22 | 23 | from ctools import obsutils |
| 23 | 24 | import sys |
| 24 | 25 | import csv |
| ... | ... | @@ -27,10 +28,10 @@ import csv |
| 27 | 28 | # ============ # |
| 28 | 29 | # cspull class # |
| 29 | 30 | # ============ # |
| 30 | -class cspull(gammalib.GApplication): | |
| 31 | +class cspull(ctools.cscript): | |
| 31 | 32 | """ |
| 32 | 33 | This class implements the pull distribution generation script. It derives |
| 33 | - from the GammaLib::GApplication class which provides support for parameter | |
| 34 | + from the ctool::cscript class which provides support for parameter | |
| 34 | 35 | files, command line arguments, and logging. In that way the Python |
| 35 | 36 | script behaves just as a regular ctool. |
| 36 | 37 | """ |
| ... | ... | @@ -52,9 +53,9 @@ class cspull(gammalib.GApplication): |
| 52 | 53 | |
| 53 | 54 | # Initialise application |
| 54 | 55 | if len(argv) == 0: |
| 55 | - gammalib.GApplication.__init__(self, self.name, self.version) | |
| 56 | + ctools.cscript.__init__(self, self.name, self.version) | |
| 56 | 57 | elif len(argv) ==1: |
| 57 | - gammalib.GApplication.__init__(self, self.name, self.version, *argv) | |
| 58 | + ctools.cscript.__init__(self, self.name, self.version, *argv) | |
| 58 | 59 | else: |
| 59 | 60 | raise TypeError("Invalid number of arguments given.") |
| 60 | 61 | |
| ... | ... | @@ -92,24 +93,28 @@ class cspull(gammalib.GApplication): |
| 92 | 93 | |
| 93 | 94 | # Create default parfile |
| 94 | 95 | pars = gammalib.GApplicationPars() |
| 96 | + pars.append(gammalib.GApplicationPar("inobs","f","h","NONE","","","Event list, counts cube, or observation definition file")) | |
| 95 | 97 | pars.append(gammalib.GApplicationPar("inmodel","f","a","$CTOOLS/share/models/crab.xml","","","Source model")) |
| 96 | 98 | pars.append(gammalib.GApplicationPar("outfile","f","a","pull.dat","","","Output file name")) |
| 97 | - pars.append(gammalib.GApplicationPar("ntrials","i","a","10","","","Number of trials")) | |
| 99 | + pars.append(gammalib.GApplicationPar("expcube","s","a","NONE","","","Exposure cube file (only needed for stacked analysis)")) | |
| 100 | + pars.append(gammalib.GApplicationPar("psfcube","s","a","NONE","","","PSF cube file (only needed for stacked analysis)")) | |
| 98 | 101 | pars.append(gammalib.GApplicationPar("caldb","s","a","dummy","","","Calibration database")) |
| 99 | 102 | pars.append(gammalib.GApplicationPar("irf","s","a","cta_dummy_irf","","","Instrument response function")) |
| 100 | - pars.append(gammalib.GApplicationPar("edisp","b","h","no","","","Apply energy dispersion?")) | |
| 101 | 103 | pars.append(gammalib.GApplicationPar("ra","r","a","83.6331","0","360","RA of pointing (deg)")) |
| 102 | 104 | pars.append(gammalib.GApplicationPar("dec","r","a","22.0145","-90","90","Dec of pointing (deg)")) |
| 103 | - pars.append(gammalib.GApplicationPar("emin","r","a","0.1","0.0","","Lower energy limit (TeV)")) | |
| 104 | - pars.append(gammalib.GApplicationPar("emax","r","a","100.0","0.0","","Upper energy limit (TeV)")) | |
| 105 | + pars.append(gammalib.GApplicationPar("emin","r","a","0.1","","","Lower energy limit (TeV)")) | |
| 106 | + pars.append(gammalib.GApplicationPar("emax","r","a","100.0","","","Upper energy limit (TeV)")) | |
| 107 | + pars.append(gammalib.GApplicationPar("tmin","r","h","0.0","","","Start time (MET in s)")) | |
| 108 | + pars.append(gammalib.GApplicationPar("tmax","r","a","1800.0","","","Duration (in s)")) | |
| 105 | 109 | pars.append(gammalib.GApplicationPar("enumbins","i","a","0","","","Number of energy bins (0=unbinned)")) |
| 106 | - pars.append(gammalib.GApplicationPar("duration","r","a","1800.0","","","Effective exposure time (s)")) | |
| 110 | + pars.append(gammalib.GApplicationPar("npix","i","h","200","","","Number of pixels for binned")) | |
| 111 | + pars.append(gammalib.GApplicationPar("binsz","r","h","0.05","","","Pixel size for binned (deg/pixel)")) | |
| 107 | 112 | pars.append(gammalib.GApplicationPar("deadc","r","h","0.95","","","Deadtime correction factor")) |
| 108 | - pars.append(gammalib.GApplicationPar("rad","r","h","5.0","","","Radius of ROI (deg)")) | |
| 113 | + pars.append(gammalib.GApplicationPar("rad","r","h","5.0","","","Radius of ROI (deg)")) | |
| 109 | 114 | pars.append(gammalib.GApplicationPar("pattern","s","h","single","","","Observation pattern (single/four)")) |
| 110 | 115 | pars.append(gammalib.GApplicationPar("offset","r","h","1.5","","","Observation pattern offset (deg)")) |
| 111 | - pars.append(gammalib.GApplicationPar("npix","i","h","200","","","Number of pixels for binned")) | |
| 112 | - pars.append(gammalib.GApplicationPar("binsz","r","h","0.05","","","Pixel size for binned (deg/pixel)")) | |
| 116 | + pars.append(gammalib.GApplicationPar("ntrials","i","a","10","","","Number of trials")) | |
| 117 | + pars.append(gammalib.GApplicationPar("edisp","b","h","no","","","Apply energy dispersion?")) | |
| 113 | 118 | pars.append_standard() |
| 114 | 119 | pars.save(parfile) |
| 115 | 120 | |
| ... | ... | @@ -121,39 +126,42 @@ class cspull(gammalib.GApplication): |
| 121 | 126 | Get parameters from parfile and setup the observation. |
| 122 | 127 | """ |
| 123 | 128 | # Get parameters |
| 124 | - if self.model == None: | |
| 125 | - self.m_inmodel = self["inmodel"].filename() | |
| 126 | - self.m_outfile = self["outfile"].filename() | |
| 127 | - self.m_ntrials = self["ntrials"].integer() | |
| 128 | - self.m_caldb = self["caldb"].string() | |
| 129 | - self.m_irf = self["irf"].string() | |
| 130 | - self.m_edisp = self["edisp"].boolean() | |
| 131 | - self.m_ra = self["ra"].real() | |
| 132 | - self.m_dec = self["dec"].real() | |
| 133 | - self.m_emin = self["emin"].real() | |
| 134 | - self.m_emax = self["emax"].real() | |
| 129 | + | |
| 130 | + # Set observation if not done before | |
| 131 | + if self.obs == None or self.obs.size() == 0: | |
| 132 | + self.obs = self.get_observations() | |
| 133 | + | |
| 134 | + # Check for requested pattern and use above | |
| 135 | + # observation parameters to set wobble pattern | |
| 136 | + self.m_pattern = self["pattern"].string() | |
| 137 | + if self.m_pattern == "four": | |
| 138 | + self.obs = self.set_obs() | |
| 139 | + | |
| 140 | + # Get number of energy bins | |
| 135 | 141 | self.m_enumbins = self["enumbins"].integer() |
| 136 | - self.m_duration = self["duration"].real() | |
| 137 | - self.m_deadc = self["deadc"].real() | |
| 138 | - self.m_rad = self["rad"].real() | |
| 139 | - self.m_pattern = self["pattern"].string() | |
| 140 | - self.m_offset = self["offset"].real() | |
| 141 | - self.m_npix = self["npix"].integer() | |
| 142 | - self.m_binsz = self["binsz"].real() | |
| 142 | + | |
| 143 | + # Read parameters for binned if requested | |
| 144 | + if not self.m_enumbins == 0: | |
| 145 | + self.m_npix = self["npix"].integer() | |
| 146 | + self.m_binsz = self["binsz"].real() | |
| 147 | + else: | |
| 148 | + # Set dummy values (required by obsutils) | |
| 149 | + self.m_npix = 0 | |
| 150 | + self.m_binsz = 0.0 | |
| 151 | + | |
| 152 | + # Set models if we have none | |
| 153 | + if self.obs.models().size() == 0: | |
| 154 | + self.obs.models(self["inmodel"].filename()) | |
| 155 | + | |
| 156 | + # Read other parameters | |
| 157 | + self.m_outfile = self["outfile"].filename() | |
| 158 | + self.m_ntrials = self["ntrials"].integer() | |
| 159 | + self.m_edisp = self["edisp"].boolean() | |
| 160 | + self.m_offset = self["offset"].real() | |
| 143 | 161 | |
| 144 | 162 | # Set some fixed parameters |
| 145 | 163 | self.m_log = False # Logging in client tools |
| 146 | 164 | self.m_debug = True # Debugging in client tools |
| 147 | - | |
| 148 | - # Setup observations | |
| 149 | - self.obs = self.set_obs() | |
| 150 | - | |
| 151 | - # Load source model | |
| 152 | - if self.m_inmodel != None: | |
| 153 | - self.model = gammalib.GModels(self.m_inmodel) | |
| 154 | - | |
| 155 | - # Append source model to observation | |
| 156 | - self.obs.models(self.model) | |
| 157 | 165 | |
| 158 | 166 | # Return |
| 159 | 167 | return |
| ... | ... | @@ -163,7 +171,7 @@ class cspull(gammalib.GApplication): |
| 163 | 171 | Set model. |
| 164 | 172 | """ |
| 165 | 173 | # Copy models |
| 166 | - self.model = models.copy() | |
| 174 | + self.obs.models(models.copy()) | |
| 167 | 175 | |
| 168 | 176 | # Return |
| 169 | 177 | return |
| ... | ... | @@ -232,18 +240,19 @@ class cspull(gammalib.GApplication): |
| 232 | 240 | |
| 233 | 241 | Keywords: |
| 234 | 242 | """ |
| 243 | + | |
| 235 | 244 | # Setup observation definition list |
| 236 | 245 | obsdeflist = obsutils.set_obs_patterns(self.m_pattern, \ |
| 237 | - ra=self.m_ra, dec=self.m_dec, \ | |
| 238 | - offset=self.m_offset) | |
| 246 | + ra=self["ra"].real(), dec=self["dec"].real(), \ | |
| 247 | + offset=self["offset"].real()) | |
| 239 | 248 | |
| 240 | 249 | # Create list of observations |
| 241 | 250 | obs = obsutils.set_obs_list(obsdeflist, \ |
| 242 | - tstart=0.0, duration=self.m_duration, \ | |
| 243 | - deadc=self.m_deadc, \ | |
| 244 | - emin=self.m_emin, emax=self.m_emax, \ | |
| 245 | - rad=self.m_rad, \ | |
| 246 | - irf=self.m_irf, caldb=self.m_caldb) | |
| 251 | + tstart=self["tmin"].real(), duration=self["tmax"].real()-self["tmin"].real(), \ | |
| 252 | + deadc=self["deadc"].real(), \ | |
| 253 | + emin=self["emin"].real(), emax=self["emax"].real(), \ | |
| 254 | + rad=self["rad"].real(), \ | |
| 255 | + irf=self["irf"].string(), caldb=self["caldb"].string()) | |
| 247 | 256 | |
| 248 | 257 | # Return observation container |
| 249 | 258 | return obs |
| ... | ... | @@ -319,7 +328,7 @@ class cspull(gammalib.GApplication): |
| 319 | 328 | |
| 320 | 329 | # Compute pull |
| 321 | 330 | fitted_value = par.value() |
| 322 | - real_value = self.model[i][k].value() | |
| 331 | + real_value = self.obs.models()[i][k].value() | |
| 323 | 332 | error = par.error() |
| 324 | 333 | if error != 0.0: |
| 325 | 334 | pull = (fitted_value - real_value) / error | ... | ... |
scripts/csresmap.py