# Licensed under a 3-clause BSD style license - see LICENSE.rst"""HESS Galactic plane survey (HGPS) catalog."""importnumpyasnpimportastropy.unitsasufromastropy.coordinatesimportAnglefromastropy.modeling.modelsimportGaussian1Dfromastropy.tableimportTablefromgammapy.estimatorsimportFluxPointsfromgammapy.mapsimportMapAxis,RegionGeomfromgammapy.modeling.modelsimportModel,Models,SkyModelfromgammapy.utils.interpolationimportScaledRegularGridInterpolatorfromgammapy.utils.scriptsimportmake_pathfromgammapy.utils.tableimporttable_row_to_dictfrom.coreimportSourceCatalog,SourceCatalogObject,format_flux_points_table__all__=["SourceCatalogHGPS","SourceCatalogLargeScaleHGPS","SourceCatalogObjectHGPS","SourceCatalogObjectHGPSComponent",]# Flux factor, used for printingFF=1e-12# Multiplicative factor to go from cm^-2 s^-1 to % Crab for integral flux > 1 TeV# Here we use the same Crab reference that's used in the HGPS paper# CRAB = crab_integral_flux(energy_min=1, reference='hess_ecpl')FLUX_TO_CRAB=100/2.26e-11FLUX_TO_CRAB_DIFF=100/3.5060459323111307e-11
[docs]classSourceCatalogObjectHGPSComponent(SourceCatalogObject):"""One Gaussian component from the HGPS catalog. See Also -------- SourceCatalogHGPS, SourceCatalogObjectHGPS """_source_name_key="Component_ID"def__init__(self,data):self.data=datadef__repr__(self):returnf"{self.__class__.__name__}({self.name!r})"def__str__(self):"""Pretty-print source data"""d=self.datass="Component {}:\n".format(d["Component_ID"])fmt="{:<20s} : {:8.3f} +/- {:.3f} deg\n"ss+=fmt.format("GLON",d["GLON"].value,d["GLON_Err"].value)ss+=fmt.format("GLAT",d["GLAT"].value,d["GLAT_Err"].value)fmt="{:<20s} : {:.3f} +/- {:.3f} deg\n"ss+=fmt.format("Size",d["Size"].value,d["Size_Err"].value)val,err=d["Flux_Map"].value,d["Flux_Map_Err"].valuefmt="{:<20s} : ({:.2f} +/- {:.2f}) x 10^-12 cm^-2 s^-1 = ({:.1f} +/- {:.1f}) % Crab"ss+=fmt.format("Flux (>1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB)returnss@propertydefname(self):"""Source name (str)"""returnself.data[self._source_name_key]
[docs]defspatial_model(self):"""Component spatial model (`~gammapy.modeling.models.GaussianSpatialModel`)."""d=self.datatag="GaussianSpatialModel"pars={"lon_0":d["GLON"],"lat_0":d["GLAT"],"sigma":d["Size"],"frame":"galactic",}errs={"lon_0":d["GLON_Err"],"lat_0":d["GLAT_Err"],"sigma":d["Size_Err"]}model=Model.create(tag,"spatial",**pars)forname,valueinerrs.items():model.parameters[name].error=valuereturnmodel
[docs]classSourceCatalogObjectHGPS(SourceCatalogObject):"""One object from the HGPS catalog. The catalog is represented by `SourceCatalogHGPS`. Examples are given there. """def__repr__(self):returnf"{self.__class__.__name__}({self.name!r})"def__str__(self):returnself.info()@propertydefflux_points(self):"""Flux points (`~gammapy.estimators.FluxPoints`)."""reference_model=SkyModel(spectral_model=self.spectral_model(),name=self.name)returnFluxPoints.from_table(self.flux_points_table,reference_model=reference_model,)
[docs]definfo(self,info="all"):"""Info string. Parameters ---------- info : {'all', 'basic', 'map', 'spec', 'flux_points', 'components', 'associations', 'id'} Comma separated list of options """ifinfo=="all":info="basic,associations,id,map,spec,flux_points,components"ss=""ops=info.split(",")if"basic"inops:ss+=self._info_basic()if"map"inops:ss+=self._info_map()if"spec"inops:ss+=self._info_spec()if"flux_points"inops:ss+=self._info_flux_points()if"components"inopsandhasattr(self,"components"):ss+=self._info_components()if"associations"inops:ss+=self._info_associations()if"id"inopsandhasattr(self,"identification_info"):ss+=self._info_id()returnss
def_info_basic(self):"""Print basic info."""d=self.datass="\n*** Basic info ***\n\n"ss+="Catalog row index (zero-based) : {}\n".format(self.row_index)ss+="{:<20s} : {}\n".format("Source name",self.name)ss+="{:<20s} : {}\n".format("Analysis reference",d["Analysis_Reference"])ss+="{:<20s} : {}\n".format("Source class",d["Source_Class"])ss+="{:<20s} : {}\n".format("Identified object",d["Identified_Object"])ss+="{:<20s} : {}\n".format("Gamma-Cat id",d["Gamma_Cat_Source_ID"])ss+="\n"returnssdef_info_associations(self):ss="\n*** Source associations info ***\n\n"lines=self.associations.pformat(max_width=-1,max_lines=-1)ss+="\n".join(lines)returnss+"\n"def_info_id(self):ss="\n*** Source identification info ***\n\n"ss+="\n".join(f"{k}: {v}"fork,vinself.identification_info.items())returnss+"\n"def_info_map(self):"""Print info from map analysis."""d=self.datass="\n*** Info from map analysis ***\n\n"ra_str=Angle(d["RAJ2000"],"deg").to_string(unit="hour",precision=0)dec_str=Angle(d["DEJ2000"],"deg").to_string(unit="deg",precision=0)ss+="{:<20s} : {:8.3f} = {}\n".format("RA",d["RAJ2000"],ra_str)ss+="{:<20s} : {:8.3f} = {}\n".format("DEC",d["DEJ2000"],dec_str)ss+="{:<20s} : {:8.3f} +/- {:.3f} deg\n".format("GLON",d["GLON"].value,d["GLON_Err"].value)ss+="{:<20s} : {:8.3f} +/- {:.3f} deg\n".format("GLAT",d["GLAT"].value,d["GLAT_Err"].value)ss+="{:<20s} : {:.3f}\n".format("Position Error (68%)",d["Pos_Err_68"])ss+="{:<20s} : {:.3f}\n".format("Position Error (95%)",d["Pos_Err_95"])ss+="{:<20s} : {:.0f}\n".format("ROI number",d["ROI_Number"])ss+="{:<20s} : {}\n".format("Spatial model",d["Spatial_Model"])ss+="{:<20s} : {}\n".format("Spatial components",d["Components"])ss+="{:<20s} : {:.1f}\n".format("TS",d["Sqrt_TS"]**2)ss+="{:<20s} : {:.1f}\n".format("sqrt(TS)",d["Sqrt_TS"])ss+="{:<20s} : {:.3f} +/- {:.3f} (UL: {:.3f}) deg\n".format("Size",d["Size"].value,d["Size_Err"].value,d["Size_UL"].value)ss+="{:<20s} : {:.3f}\n".format("R70",d["R70"])ss+="{:<20s} : {:.3f}\n".format("RSpec",d["RSpec"])ss+="{:<20s} : {:.1f}\n".format("Total model excess",d["Excess_Model_Total"])ss+="{:<20s} : {:.1f}\n".format("Excess in RSpec",d["Excess_RSpec"])ss+="{:<20s} : {:.1f}\n".format("Model Excess in RSpec",d["Excess_RSpec_Model"])ss+="{:<20s} : {:.1f}\n".format("Background in RSpec",d["Background_RSpec"])ss+="{:<20s} : {:.1f} hours\n".format("Livetime",d["Livetime"].value)ss+="{:<20s} : {:.2f}\n".format("Energy threshold",d["Energy_Threshold"])val,err=d["Flux_Map"].value,d["Flux_Map_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: 501"Source flux (>1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB,)ss+="\nFluxes in RSpec (> 1 TeV):\n"ss+="{:<30s} : {:.3f} x 10^-12 cm^-2 s^-1 = {:5.2f} % Crab\n".format("Map measurement",d["Flux_Map_RSpec_Data"].value/FF,d["Flux_Map_RSpec_Data"].value*FLUX_TO_CRAB,)ss+="{:<30s} : {:.3f} x 10^-12 cm^-2 s^-1 = {:5.2f} % Crab\n".format("Source model",d["Flux_Map_RSpec_Source"].value/FF,d["Flux_Map_RSpec_Source"].value*FLUX_TO_CRAB,)ss+="{:<30s} : {:.3f} x 10^-12 cm^-2 s^-1 = {:5.2f} % Crab\n".format("Other component model",d["Flux_Map_RSpec_Other"].value/FF,d["Flux_Map_RSpec_Other"].value*FLUX_TO_CRAB,)ss+="{:<30s} : {:.3f} x 10^-12 cm^-2 s^-1 = {:5.2f} % Crab\n".format("Large scale component model",d["Flux_Map_RSpec_LS"].value/FF,d["Flux_Map_RSpec_LS"].value*FLUX_TO_CRAB,)ss+="{:<30s} : {:.3f} x 10^-12 cm^-2 s^-1 = {:5.2f} % Crab\n".format("Total model",d["Flux_Map_RSpec_Total"].value/FF,d["Flux_Map_RSpec_Total"].value*FLUX_TO_CRAB,)ss+="{:<35s} : {:5.1f} %\n".format("Containment in RSpec",100*d["Containment_RSpec"])ss+="{:<35s} : {:5.1f} %\n".format("Contamination in RSpec",100*d["Contamination_RSpec"])label,val=("Flux correction (RSpec -> Total)",100*d["Flux_Correction_RSpec_To_Total"],)ss+=f"{label:<35s} : {val:5.1f} %\n"label,val=("Flux correction (Total -> RSpec)",100*(1/d["Flux_Correction_RSpec_To_Total"]),)ss+=f"{label:<35s} : {val:5.1f} %\n"returnssdef_info_spec(self):"""Print info from spectral analysis."""d=self.datass="\n*** Info from spectral analysis ***\n\n"ss+="{:<20s} : {:.1f} hours\n".format("Livetime",d["Livetime_Spec"].value)lo=d["Energy_Range_Spec_Min"].valuehi=d["Energy_Range_Spec_Max"].valuess+="{:<20s} : {:.2f} to {:.2f} TeV\n".format("Energy range:",lo,hi)ss+="{:<20s} : {:.1f}\n".format("Background",d["Background_Spec"])ss+="{:<20s} : {:.1f}\n".format("Excess",d["Excess_Spec"])ss+="{:<20s} : {}\n".format("Spectral model",d["Spectral_Model"])val=d["TS_ECPL_over_PL"]ss+="{:<20s} : {:.1f}\n".format("TS ECPL over PL",val)val=d["Flux_Spec_Int_1TeV"].valueerr=d["Flux_Spec_Int_1TeV_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: E501"Best-fit model flux(> 1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB,)val=d["Flux_Spec_Energy_1_10_TeV"].valueerr=d["Flux_Spec_Energy_1_10_TeV_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 erg cm^-2 s^-1\n".format("Best-fit model energy flux(1 to 10 TeV)",val/FF,err/FF)ss+=self._info_spec_pl()ss+=self._info_spec_ecpl()returnssdef_info_spec_pl(self):d=self.datass="{:<20s} : {:.2f}\n".format("Pivot energy",d["Energy_Spec_PL_Pivot"])val=d["Flux_Spec_PL_Diff_Pivot"].valueerr=d["Flux_Spec_PL_Diff_Pivot_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 TeV^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: E501"Flux at pivot energy",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB_DIFF,)val=d["Flux_Spec_PL_Int_1TeV"].valueerr=d["Flux_Spec_PL_Int_1TeV_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: E501"PL Flux(> 1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB,)val=d["Flux_Spec_PL_Diff_1TeV"].valueerr=d["Flux_Spec_PL_Diff_1TeV_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 TeV^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: E501"PL Flux(@ 1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB_DIFF,)val=d["Index_Spec_PL"]err=d["Index_Spec_PL_Err"]ss+="{:<20s} : {:.2f} +/- {:.2f}\n".format("PL Index",val,err)returnssdef_info_spec_ecpl(self):d=self.datass=""# ss = '{:<20s} : {:.1f}\n'.format('Pivot energy', d['Energy_Spec_ECPL_Pivot'])val=d["Flux_Spec_ECPL_Diff_1TeV"].valueerr=d["Flux_Spec_ECPL_Diff_1TeV_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 TeV^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: E501"ECPL Flux(@ 1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB_DIFF,)val=d["Flux_Spec_ECPL_Int_1TeV"].valueerr=d["Flux_Spec_ECPL_Int_1TeV_Err"].valuess+="{:<20s} : ({:.3f} +/- {:.3f}) x 10^-12 cm^-2 s^-1 = ({:.2f} +/- {:.2f}) % Crab\n".format(# noqa: E501"ECPL Flux(> 1 TeV)",val/FF,err/FF,val*FLUX_TO_CRAB,err*FLUX_TO_CRAB,)val=d["Index_Spec_ECPL"]err=d["Index_Spec_ECPL_Err"]ss+="{:<20s} : {:.2f} +/- {:.2f}\n".format("ECPL Index",val,err)val=d["Lambda_Spec_ECPL"].valueerr=d["Lambda_Spec_ECPL_Err"].valuess+="{:<20s} : {:.3f} +/- {:.3f} TeV^-1\n".format("ECPL Lambda",val,err)# Use Gaussian analytical error propagation,# tested against the uncertainties packageerr=err/val**2val=1.0/valss+="{:<20s} : {:.2f} +/- {:.2f} TeV\n".format("ECPL E_cut",val,err)returnssdef_info_flux_points(self):"""Print flux point results"""d=self.datass="\n*** Flux points info ***\n\n"ss+="Number of flux points: {}\n".format(d["N_Flux_Points"])ss+="Flux points table: \n\n"lines=format_flux_points_table(self.flux_points_table).pformat(max_width=-1,max_lines=-1)ss+="\n".join(lines)returnss+"\n"def_info_components(self):"""Print info about the components."""ss="\n*** Gaussian component info ***\n\n"ss+="Number of components: {}\n".format(len(self.components))ss+="{:<20s} : {}\n\n".format("Spatial components",self.data["Components"])forcomponentinself.components:ss+=str(component)ss+="\n\n"returnss@propertydefenergy_range(self):"""Spectral model energy range (`~astropy.units.Quantity` with length 2)."""energy_min,energy_max=(self.data["Energy_Range_Spec_Min"],self.data["Energy_Range_Spec_Max"],)ifnp.isnan(energy_min):energy_min=u.Quantity(0.2,"TeV")ifnp.isnan(energy_max):energy_max=u.Quantity(50,"TeV")returnu.Quantity([energy_min,energy_max],"TeV")
[docs]defspectral_model(self,which="best"):"""Spectral model (`~gammapy.modeling.models.SpectralModel`). One of the following models (given by ``Spectral_Model`` in the catalog): - ``PL`` : `~gammapy.modeling.models.PowerLawSpectralModel` - ``ECPL`` : `~gammapy.modeling.models.ExpCutoffPowerLawSpectralModel` Parameters ---------- which : {'best', 'pl', 'ecpl'} Which spectral model """data=self.dataifwhich=="best":spec_type=self.data["Spectral_Model"].strip().lower()elifwhichin{"pl","ecpl"}:spec_type=whichelse:raiseValueError(f"Invalid selection: which = {which!r}")ifspec_type=="pl":tag="PowerLawSpectralModel"pars={"index":data["Index_Spec_PL"],"amplitude":data["Flux_Spec_PL_Diff_Pivot"],"reference":data["Energy_Spec_PL_Pivot"],}errs={"amplitude":data["Flux_Spec_PL_Diff_Pivot_Err"],"index":data["Index_Spec_PL_Err"],}elifspec_type=="ecpl":tag="ExpCutoffPowerLawSpectralModel"pars={"index":data["Index_Spec_ECPL"],"amplitude":data["Flux_Spec_ECPL_Diff_Pivot"],"reference":data["Energy_Spec_ECPL_Pivot"],"lambda_":data["Lambda_Spec_ECPL"],}errs={"index":data["Index_Spec_ECPL_Err"],"amplitude":data["Flux_Spec_ECPL_Diff_Pivot_Err"],"lambda_":data["Lambda_Spec_ECPL_Err"],}else:raiseValueError(f"Invalid spec_type: {spec_type}")model=Model.create(tag,"spectral",**pars)errs["reference"]=0*u.TeVforname,valueinerrs.items():model.parameters[name].error=valuereturnmodel
[docs]defspatial_model(self):"""Spatial model (`~gammapy.modeling.models.SpatialModel`). One of the following models (given by ``Spatial_Model`` in the catalog): - ``Point-Like`` or has a size upper limit : `~gammapy.modeling.models.PointSpatialModel` - ``Gaussian``: `~gammapy.modeling.models.GaussianSpatialModel` - ``2-Gaussian`` or ``3-Gaussian``: composite model (using ``+`` with Gaussians) - ``Shell``: `~gammapy.modeling.models.ShellSpatialModel` """d=self.datapars={"lon_0":d["GLON"],"lat_0":d["GLAT"],"frame":"galactic"}errs={"lon_0":d["GLON_Err"],"lat_0":d["GLAT_Err"]}spatial_type=self.data["Spatial_Model"]ifspatial_type=="Point-Like":tag="PointSpatialModel"elifspatial_type=="Gaussian":tag="GaussianSpatialModel"pars["sigma"]=d["Size"]errs["sigma"]=d["Size_Err"]elifspatial_typein{"2-Gaussian","3-Gaussian"}:raiseValueError("For Gaussian or Multi-Gaussian models, use sky_model()!")elifspatial_type=="Shell":# HGPS contains no information on shell width# Here we assume a 5% shell width for all shells.tag="ShellSpatialModel"pars["radius"]=0.95*d["Size"]pars["width"]=d["Size"]-pars["radius"]errs["radius"]=d["Size_Err"]else:raiseValueError(f"Invalid spatial_type: {spatial_type}")model=Model.create(tag,"spatial",**pars)forname,valueinerrs.items():model.parameters[name].error=valuereturnmodel
[docs]defsky_model(self,which="best"):"""Source sky model. Parameters ---------- which : {'best', 'pl', 'ecpl'} Which spectral model Returns ------- sky_model : `~gammapy.modeling.models.Models` Models of the catalog object. """models=self.components_models(which=which)iflen(models)>1:geom=self._get_components_geom(models)returnmodels.to_template_sky_model(geom=geom,name=self.name)else:returnmodels[0]
[docs]defcomponents_models(self,which="best",linked=False):"""Models of the source components. Parameters ---------- which : {'best', 'pl', 'ecpl'} Which spectral model linked : bool Each sub-component of a source is given as a different `SkyModel` If True the spectral parameters except the mormalisation are linked. Default is False Returns ------- sky_model : `~gammapy.modeling.models.Models` Models of the catalog object. """spatial_type=self.data["Spatial_Model"]missing_size=(spatial_type=="Gaussian"andself.spatial_model().sigma.value==0)ifspatial_typein{"2-Gaussian","3-Gaussian"}ormissing_size:models=[]spectral_model=self.spectral_model(which=which)forcomponentinself.components:spec_component=spectral_model.copy()weight=component.data["Flux_Map"]/self.data["Flux_Map"]spec_component.parameters["amplitude"].value*=weightiflinked:fornameinspec_component.parameters.names:ifnamenotin["norm","amplitude"]:spec_component.__dict__[name]=spectral_model.parameters[name]model=SkyModel(spatial_model=component.spatial_model(),spectral_model=spec_component,name=component.name,)models.append(model)else:models=[SkyModel(spatial_model=self.spatial_model(),spectral_model=self.spectral_model(which=which),name=self.name,)]returnModels(models)
[docs]defto_models(self,which="best",components_status="independent"):"""Create Models object from catalogue Parameters ---------- which : {'best', 'pl', 'ecpl'} Which spectral model components_status : {'independent', 'linked', 'merged'} Relation between the sources components: 'independent' : each sub-component of a source is given as a different `SkyModel` (Default) 'linked' : each sub-component of a source is given as a different `SkyModel` but the spectral parameters except the mormalisation are linked. 'merged' : the sub-components are merged into a single `SkyModel` given as a `~gammapy.modeling.models.TemplateSpatialModel` with a `~gammapy.modeling.models.PowerLawNormSpectralModel`. In that case the relavtie weights between the components cannot be adjusted. Returns ------- models : `~gammapy.modeling.models.Models` Models of the catalog. """models=[]forsourceinself:ifcomponents_status=="merged":m=[source.sky_model(which=which)]else:m=source.components_models(which=which,linked=components_status=="linked")models.extend(m)returnModels(models)
[docs]classSourceCatalogLargeScaleHGPS:"""Gaussian band model. This 2-dimensional model is Gaussian in ``y`` for a given ``x``, and the Gaussian parameters can vary in ``x``. One application of this model is the diffuse emission along the Galactic plane, i.e. ``x = GLON`` and ``y = GLAT``. Parameters ---------- table : `~astropy.table.Table` Table of Gaussian parameters. ``x``, ``amplitude``, ``mean``, ``stddev``. interp_kwargs : dict Keyword arguments passed to `ScaledRegularGridInterpolator` """def__init__(self,table,interp_kwargs=None):interp_kwargs=interp_kwargsor{}interp_kwargs.setdefault("values_scale","lin")self.table=tableglon=Angle(self.table["GLON"]).wrap_at("180d")interps={}forcolumnintable.colnames:values=self.table[column].quantityinterp=ScaledRegularGridInterpolator((glon,),values,**interp_kwargs)interps[column]=interpself._interp=interpsdef_interpolate_parameter(self,parname,glon):glon=glon.wrap_at("180d")returnself._interp[parname]((np.asanyarray(glon),),clip=False)
[docs]defpeak_brightness(self,glon):"""Peak brightness at a given longitude. Parameters ---------- glon : `~astropy.coordinates.Angle` Galactic Longitude. """returnself._interpolate_parameter("Surface_Brightness",glon)
[docs]defpeak_brightness_error(self,glon):"""Peak brightness error at a given longitude. Parameters ---------- glon : `~astropy.coordinates.Angle` Galactic Longitude. """returnself._interpolate_parameter("Surface_Brightness_Err",glon)
[docs]defwidth(self,glon):"""Width at a given longitude. Parameters ---------- glon : `~astropy.coordinates.Angle` Galactic Longitude. """returnself._interpolate_parameter("Width",glon)
[docs]defwidth_error(self,glon):"""Width error at a given longitude. Parameters ---------- glon : `~astropy.coordinates.Angle` Galactic Longitude. """returnself._interpolate_parameter("Width_Err",glon)
[docs]defpeak_latitude(self,glon):"""Peak position at a given longitude. Parameters ---------- glon : `~astropy.coordinates.Angle` Galactic Longitude. """returnself._interpolate_parameter("GLAT",glon)
[docs]defpeak_latitude_error(self,glon):"""Peak position error at a given longitude. Parameters ---------- glon : `~astropy.coordinates.Angle` Galactic Longitude. """returnself._interpolate_parameter("GLAT_Err",glon)
[docs]defevaluate(self,position):"""Evaluate model at a given position. Parameters ---------- position : `~astropy.coordinates.SkyCoord` Position on the sky. """glon,glat=position.galactic.l,position.galactic.bwidth=self.width(glon)amplitude=self.peak_brightness(glon)mean=self.peak_latitude(glon)returnGaussian1D.evaluate(glat,amplitude=amplitude,mean=mean,stddev=width)