lsst.sims.maf.slicers package¶
Submodules¶
lsst.sims.maf.slicers.baseSlicer module¶
-
class
lsst.sims.maf.slicers.baseSlicer.
BaseSlicer
(verbose=True, badval=- 666)[source]¶ Bases:
object
Base class for all slicers: sets required methods and implements common functionality.
After first construction, the slicer should be ready for setupSlicer to define slicePoints, which will let the slicer ‘slice’ data and generate plots. After init after a restore: everything necessary for using slicer for plotting or saving/restoring metric data should be present (although slicer does not need to be able to slice data again and generally will not be able to).
- Parameters
-
__next__
()[source]¶ Returns results of self._sliceSimData when iterating over slicer.
Results of self._sliceSimData should be dictionary of {‘idxs’: the data indexes relevant for this slice of the slicer, ‘slicePoint’: the metadata for the slicePoint, which always includes ‘sid’ key for ID of slicePoint.}
-
outputJSON
(metricValues, metricName='', simDataName='', metadata='', plotDict=None)[source]¶ Send metric data to JSON streaming API, along with a little bit of metadata.
This method will only work for metrics where the metricDtype is float or int, as JSON will not interpret more complex data properly. These values can’t be plotted anyway though.
- Parameters
metricValues (np.ma.MaskedArray or np.ndarray) – The metric values.
metricName (str, optional) – The name of the metric. Default ‘’.
simDataName (str, optional) – The name of the simulated data source. Default ‘’.
metadata (str, optional) – The metadata about this metric. Default ‘’.
plotDict (dict, optional.) – The plotDict for this metric bundle. Default None.
- Returns
StringIO object containing a header dictionary with metricName/metadata/simDataName/slicerName, and plot labels from plotDict, and metric values/data for plot. if oneDSlicer, the data is [ [bin_left_edge, value], [bin_left_edge, value]..]. if a spatial slicer, the data is [ [lon, lat, value], [lon, lat, value] ..].
- Return type
StringIO
-
readData
(infilename)[source]¶ Read metric data from disk, along with the info to rebuild the slicer (minus new slicing capability).
-
registry
¶
-
setupSlicer
(simData, maps=None)[source]¶ Set up Slicer for data slicing.
Set up internal parameters necessary for slicer to slice data and generates indexes on simData. Also sets _sliceSimData for a particular slicer.
- Parameters
simData (np.recarray) – The simulated data to be sliced.
maps (list of lsst.sims.maf.maps objects, optional.) – Maps to apply at each slicePoint, to add to the slicePoint metadata. Default None.
-
writeData
(outfilename, metricValues, metricName='', simDataName='', constraint=None, metadata='', plotDict=None, displayDict=None)[source]¶ Save metric values along with the information required to re-build the slicer.
- Parameters
outfilename (str) – The output file name.
metricValues (np.ma.MaskedArray or np.ndarray) – The metric values to save to disk.
lsst.sims.maf.slicers.baseSpatialSlicer module¶
-
class
lsst.sims.maf.slicers.baseSpatialSlicer.
BaseSpatialSlicer
(lonCol='fieldRA', latCol='fieldDec', latLonDeg=True, verbose=True, badval=- 666, leafsize=100, radius=1.75, useCamera=False, rotSkyPosColName='rotSkyPos', mjdColName='observationStartMJD', chipNames='all', scienceChips=True)[source]¶ Bases:
lsst.sims.maf.slicers.baseSlicer.BaseSlicer
Base spatial slicer object, contains additional functionality for spatial slicing, including setting up and traversing a kdtree containing the simulated data points.
- Parameters
lonCol (str, optional) – Name of the longitude (RA equivalent) column to use from the input data. Default fieldRA
latCol (str, optional) – Name of the latitude (Dec equivalent) column to use from the input data. Default fieldDec
latLonDeg (boolean, optional) – Flag indicating whether lat and lon values from input data are in degrees (True) or radians (False). Default True.
verbose (boolean, optional) – Flag to indicate whether or not to write additional information to stdout during runtime. Default True.
badval (float, optional) – Bad value flag, relevant for plotting. Default -666.
leafsize (int, optional) – Leafsize value for kdtree. Default 100.
radius (float, optional) – Radius for matching in the kdtree. Equivalent to the radius of the FOV. Degrees. Default 1.75.
useCamera (boolean, optional) – Flag to indicate whether to use the LSST camera footprint or not. Default False.
rotSkyPosColName (str, optional) – Name of the rotSkyPos column in the input data. Only used if useCamera is True. Describes the orientation of the camera orientation compared to the sky. Default rotSkyPos.
mjdColName (str, optional) – Name of the exposure time column. Only used if useCamera is True. Default observationStartMJD.
chipNames (array-like, optional) – List of chips to accept, if useCamera is True. This lets users turn ‘on’ only a subset of chips. Default ‘all’ - this uses all chips in the camera.
scienceChips (bool (True)) – Do not include wavefront sensors when checking if a point landed on a chip.
-
setupSlicer
(simData, maps=None)[source]¶ Use simData[self.lonCol] and simData[self.latCol] (in radians) to set up KDTree.
- Parameters
simData (numpy.recarray) – The simulated data, including the location of each pointing.
maps (list of lsst.sims.maf.maps objects, optional) – List of maps (such as dust extinction) that will run to build up additional metadata at each slicePoint. This additional metadata is available to metrics via the slicePoint dictionary. Default None.
lsst.sims.maf.slicers.healpixComCamSlicer module¶
-
class
lsst.sims.maf.slicers.healpixComCamSlicer.
HealpixComCamSlicer
(nside=128, lonCol='fieldRA', latCol='fieldDec', latLonDeg=True, verbose=True, badval=- 1.6375e+30, useCache=True, leafsize=100, useCamera=False, rotSkyPosColName='rotSkyPos', mjdColName='observationStartMJD', chipNames=['R:2,2 S:0,0', 'R:2,2 S:0,1', 'R:2,2 S:0,2', 'R:2,2 S:1,0', 'R:2,2 S:1,1', 'R:2,2 S:1,2', 'R:2,2 S:2,0', 'R:2,2 S:2,1', 'R:2,2 S:2,2'], side_length=0.7)[source]¶ Bases:
lsst.sims.maf.slicers.healpixSlicer.HealpixSlicer
Slicer that uses the ComCam footprint to decide if observations overlap a healpixel center
-
setupSlicer
(simData, maps=None)[source]¶ Use simData[self.lonCol] and simData[self.latCol] (in radians) to set up KDTree.
- Parameters
simData (numpy.recarray) – The simulated data, including the location of each pointing.
maps (list of lsst.sims.maf.maps objects, optional) – List of maps (such as dust extinction) that will run to build up additional metadata at each slicePoint. This additional metadata is available to metrics via the slicePoint dictionary. Default None.
-
lsst.sims.maf.slicers.healpixSDSSSlicer module¶
-
class
lsst.sims.maf.slicers.healpixSDSSSlicer.
HealpixSDSSSlicer
(nside=128, lonCol='RA1', latCol='Dec1', verbose=True, useCache=True, radius=0.2833333333333333, leafsize=100, **kwargs)[source]¶ Bases:
lsst.sims.maf.slicers.healpixSlicer.HealpixSlicer
For use with SDSS stripe 82 square images
lsst.sims.maf.slicers.healpixSlicer module¶
A slicer that uses a Healpix grid to calculate metric values (at the center of each healpixel).
-
class
lsst.sims.maf.slicers.healpixSlicer.
HealpixSlicer
(nside=128, lonCol='fieldRA', latCol='fieldDec', latLonDeg=True, verbose=True, badval=- 1.6375e+30, useCache=True, leafsize=100, radius=1.75, useCamera=False, rotSkyPosColName='rotSkyPos', mjdColName='observationStartMJD', chipNames='all')[source]¶ Bases:
lsst.sims.maf.slicers.baseSpatialSlicer.BaseSpatialSlicer
A spatial slicer that evaluates pointings on a healpix-based grid.
Note that a healpix slicer is intended to evaluate the sky on a spatial grid. Usually this grid will be something like RA as the lonCol and Dec as the latCol. However, it could also be altitude and azimuth, in which case altitude as latCol, and azimuth as lonCol. An additional alternative is to use HA/Dec, with lonCol=HA, latCol=Dec.
When plotting with RA/Dec, the default HealpixSkyMap can be used, corresponding to {‘rot’: (0, 0, 0), ‘flip’: ‘astro’}. When plotting with alt/az, either the LambertSkyMap can be used (if Basemap is installed) or the HealpixSkyMap can be used with a modified plotDict, {‘rot’: (90, 90, 90), ‘flip’: ‘geo’}. When plotting with HA/Dec, only the HealpixSkyMap can be used, with a modified plotDict of {‘rot’: (0, -30, 0), ‘flip’: ‘geo’}.
- Parameters
nside (int, optional) – The nside parameter of the healpix grid. Must be a power of 2. Default 128.
lonCol (str, optional) – Name of the longitude (RA equivalent) column to use from the input data. Default fieldRA
latCol (str, optional) – Name of the latitude (Dec equivalent) column to use from the input data. Default fieldDec
latLonDeg (boolean, optional) – Flag indicating whether the lat and lon values in the input data are in degrees (True) or radians (False). Default True.
verbose (boolean, optional) – Flag to indicate whether or not to write additional information to stdout during runtime. Default True.
badval (float, optional) – Bad value flag, relevant for plotting. Default the hp.UNSEEN value (in order to properly flag bad data points for plotting with the healpix plotting routines). This should not be changed.
useCache (boolean) – Flag allowing the user to indicate whether or not to cache (and reuse) metric results calculated with the same set of simulated data pointings. This can be safely set to True for slicers not using maps and will result in increased speed. When calculating metric results using maps, the metadata at each individual ra/dec point may influence the metric results and so useCache should be set to False. Default True.
leafsize (int, optional) – Leafsize value for kdtree. Default 100.
radius (float, optional) – Radius for matching in the kdtree. Equivalent to the radius of the FOV. Degrees. Default 1.75.
useCamera (boolean, optional) – Flag to indicate whether to use the LSST camera footprint or not. Default False.
rotSkyPosColName (str, optional) – Name of the rotSkyPos column in the input data. Only used if useCamera is True. Describes the orientation of the camera orientation compared to the sky. Default rotSkyPos.
mjdColName (str, optional) – Name of the exposure time column. Only used if useCamera is True. Default observationStartMJD.
chipNames (array-like, optional) – List of chips to accept, if useCamera is True. This lets users turn ‘on’ only a subset of chips. Default ‘all’ - this uses all chips in the camera.
lsst.sims.maf.slicers.healpixSubsetSlicer module¶
A HealpixSubsetSlicer - define the subset of healpixels to use to calculate metrics.
-
class
lsst.sims.maf.slicers.healpixSubsetSlicer.
HealpixSubsetSlicer
(nside, hpid, lonCol='fieldRA', latCol='fieldDec', latLonDeg=True, verbose=True, badval=- 1.6375e+30, useCache=True, leafsize=100, radius=1.75, useCamera=False, rotSkyPosColName='rotSkyPos', mjdColName='observationStartMJD', chipNames='all')[source]¶ Bases:
lsst.sims.maf.slicers.baseSpatialSlicer.BaseSpatialSlicer
A spatial slicer that evaluates pointings on a subset of a healpix-based grid. The advantage of using this healpixSubsetSlicer (rather than just putting the RA/Dec values into the UserPointsSlicer, which is another valid approach) is that you preserve the full healpix array. This means you could do things like calculate the power spectrum and plot without remapping into healpixels first. The downside is that you must first (externally) define the healpixels that you wish to use - the lsst.sims.featureScheduler.footprints is a useful add-on here.
When plotting with RA/Dec, the default HealpixSkyMap can be used, corresponding to {‘rot’: (0, 0, 0), ‘flip’: ‘astro’}.
- Parameters
nside (int) – The nside parameter of the healpix grid. Must be a power of 2.
hpid (np.ndarray) – The subset of healpix id’s to use to calculate the metric. Because the hpid should be defined based on a particular nside, these first two arguments are not optional for this slicer.
lonCol (str, optional) – Name of the longitude (RA equivalent) column to use from the input data. Default fieldRA
latCol (str, optional) – Name of the latitude (Dec equivalent) column to use from the input data. Default fieldDec
latLonDeg (boolean, optional) – Flag indicating whether the lat and lon values in the input data are in degrees (True) or radians (False). Default True.
verbose (boolean, optional) – Flag to indicate whether or not to write additional information to stdout during runtime. Default True.
badval (float, optional) – Bad value flag, relevant for plotting. Default the hp.UNSEEN value (in order to properly flag bad data points for plotting with the healpix plotting routines). This should not be changed.
useCache (boolean) – Flag allowing the user to indicate whether or not to cache (and reuse) metric results calculated with the same set of simulated data pointings. This can be safely set to True for slicers not using maps and will result in increased speed. When calculating metric results using maps, the metadata at each individual ra/dec point may influence the metric results and so useCache should be set to False. Default True.
leafsize (int, optional) – Leafsize value for kdtree. Default 100.
radius (float, optional) – Radius for matching in the kdtree. Equivalent to the radius of the FOV. Degrees. Default 1.75.
useCamera (boolean, optional) – Flag to indicate whether to use the LSST camera footprint or not. Default False.
rotSkyPosColName (str, optional) – Name of the rotSkyPos column in the input data. Only used if useCamera is True. Describes the orientation of the camera orientation compared to the sky. Default rotSkyPos.
mjdColName (str, optional) – Name of the exposure time column. Only used if useCamera is True. Default observationStartMJD.
chipNames (array-like, optional) – List of chips to accept, if useCamera is True. This lets users turn ‘on’ only a subset of chips. Default ‘all’ - this uses all chips in the camera.
-
setupSlicer
(simData, maps=None)[source]¶ Use simData[self.lonCol] and simData[self.latCol] (in radians) to set up KDTree.
- Parameters
simData (numpy.recarray) – The simulated data, including the location of each pointing.
maps (list of lsst.sims.maf.maps objects, optional) – List of maps (such as dust extinction) that will run to build up additional metadata at each slicePoint. This additional metadata is available to metrics via the slicePoint dictionary. Default None.
lsst.sims.maf.slicers.hourglassSlicer module¶
-
class
lsst.sims.maf.slicers.hourglassSlicer.
HourglassSlicer
(verbose=True, badval=- 666)[source]¶ Bases:
lsst.sims.maf.slicers.uniSlicer.UniSlicer
Slicer to make the filter hourglass plots
lsst.sims.maf.slicers.moSlicer module¶
-
class
lsst.sims.maf.slicers.moSlicer.
MoObjSlicer
(Hrange=None, verbose=True, badval=0)[source]¶ Bases:
lsst.sims.maf.slicers.baseSlicer.BaseSlicer
Slice moving object _observations_, per object and optionally clone/per H value.
Iteration over the MoObjSlicer will go as: * iterate over each orbit; * if Hrange is not None, for each orbit, iterate over Hrange.
- Parameters
Hrange (numpy.ndarray or None) – The H values to clone the orbital parameters over. If Hrange is None, will not clone orbits.
-
readObs
(obsFile)[source]¶ Read observations of the solar system objects (such as created by sims_movingObjects).
- Parameters
obsFile (str) – The file containing the observation information.
-
setupSlicer
(orbitFile, delim=None, skiprows=None, obsFile=None)[source]¶ Set up the slicer and read orbitFile and obsFile from disk.
Sets self.orbits (with orbit parameters), self.allObs, and self.obs self.orbitFile and self.obsFile
- Parameters
orbitFile (str) – The file containing the orbit information. This is necessary, in order to be able to generate plots.
obsFile (str, optional) – The file containing the observations of each object, optional. If not provided (default, None), then the slicer will not be able to ‘slice’, but can still plot.
lsst.sims.maf.slicers.movieSlicer module¶
-
class
lsst.sims.maf.slicers.movieSlicer.
MovieSlicer
(sliceColName=None, sliceColUnits=None, bins=None, binMin=None, binMax=None, binsize=None, verbose=True, badval=0, cumulative=True, forceNoFfmpeg=False)[source]¶
lsst.sims.maf.slicers.nDSlicer module¶
lsst.sims.maf.slicers.oneDSlicer module¶
-
class
lsst.sims.maf.slicers.oneDSlicer.
OneDSlicer
(sliceColName=None, sliceColUnits=None, bins=None, binMin=None, binMax=None, binsize=None, verbose=True, badval=0)[source]¶ Bases:
lsst.sims.maf.slicers.baseSlicer.BaseSlicer
oneD Slicer.
lsst.sims.maf.slicers.opsimFieldSlicer module¶
-
class
lsst.sims.maf.slicers.opsimFieldSlicer.
OpsimFieldSlicer
(simDataFieldIdColName='fieldId', simDataFieldRaColName='fieldRA', simDataFieldDecColName='fieldDec', latLonDeg=False, fieldIdColName='fieldId', fieldRaColName='fieldRA', fieldDecColName='fieldDec', verbose=True, badval=- 666)[source]¶ Bases:
lsst.sims.maf.slicers.baseSpatialSlicer.BaseSpatialSlicer
A spatial slicer that evaluates pointings based on matched IDs between the simData and fieldData.
Note that this slicer uses the fieldId of the simulated data fields to generate the spatial matches. Thus, it is not suitable for use in evaluating dithering or high resolution metrics (use the HealpixSlicer instead for those use-cases).
When the slicer is set up, it takes two arrays: fieldData and simData. FieldData is a numpy.recarray containing the information about the fields - this is the basis for slicing. The simData is a numpy.recarray that holds the information about the pointings - this is the data that is matched against the fieldData.
- Parameters
simDataFieldIDColName (str, optional) – Name of the column in simData for the fieldId Default fieldId.
simDataFieldRaColName (str, optional) – Name of the column in simData for the RA. Default fieldRA.
simDataFieldDecColName (str, optional) – Name of the column in simData for the fieldDec. Default fieldDec.
latLongDeg (bool, optional) – Whether the RA/Dec values in fieldData are in degrees. If using a standard metricBundleGroup to run the metric, FieldData is fetched by utils.getFieldData, which always returns radians (so the default here is False).
fieldIdColName (str, optional) – Name of the column in the fieldData for the fieldId (to match with simData). Default fieldId.
fieldRaColName (str, optional) – Name of the column in the fieldData for the RA (used for plotting). Default fieldRA.
fieldDecColName (str, optional) – Name of the column in the fieldData for the Dec (used for plotting). Default fieldDec.
verbose (boolean, optional) – Flag to indicate whether or not to write additional information to stdout during runtime. Default True.
badval (float, optional) – Bad value flag, relevant for plotting. Default -666.
-
setupSlicer
(simData, fieldData, maps=None)[source]¶ Set up opsim field slicer object.
- Parameters
simData (numpy.recarray) – Contains the simulation pointing history.
fieldData (numpy.recarray) – Contains the field information (ID, Ra, Dec) about how to slice the simData. For example, only fields in the fieldData table will be matched against the simData. RA and Dec should be in degrees.
maps (list of lsst.sims.maf.maps objects, optional) – Maps to run and provide additional metadata at each slicePoint. Default None.
lsst.sims.maf.slicers.orbits module¶
-
class
lsst.sims.maf.slicers.orbits.
Orbits
[source]¶ Bases:
object
Orbits reads, checks for required values, and stores orbit parameters for moving objects.
Instantiate the class and then use readOrbits or setOrbits to set the orbit values.
self.orbits stores the orbital parameters, as a pandas dataframe. self.dataCols defines the columns required, although objId, H, g, and sed_filename are optional.
-
assignSed
(orbits, randomSeed=None)[source]¶ Assign either a C or S type SED, depending on the semi-major axis of the object. P(C type) = 0 (a<2); 0.5*a - 1 (2<a<4); 1 (a > 4), based on figure 23 from Ivezic et al 2001 (AJ, 122, 2749).
- Parameters
orbits (pandas.DataFrame, pandas.Series or numpy.ndarray) – Array-like object containing orbital parameter information.
- Returns
Array containing the SED type for each object in ‘orbits’.
- Return type
numpy.ndarray
-
readOrbits
(orbitfile, delim=None, skiprows=None)[source]¶ Read orbits from a file, generating a pandas dataframe containing columns matching dataCols, for the appropriate orbital parameter format (currently accepts COM, KEP or CAR formats).
After reading and standardizing the column names, calls self.setOrbits to validate the orbital parameters. Expects angles in orbital element formats to be in degrees.
Note that readOrbits uses pandas.read_csv to read the data file with the orbital parameters. Thus, it should have column headers specifying the column names .. unless skiprows = -1 or there is just no header line at all. in which case it is assumed to be a standard DES format file, with no header line.
- Parameters
orbitfile (str) – The name of the input file containing orbital parameter information.
delim (str, optional) – The delimiter for the input orbit file. Default is None, will use delim_whitespace=True.
skiprows (int, optional) – The number of rows to skip before reading the header information for pandas. Default is None, which will trigger a check of the file to look for the header columns.
-
setOrbits
(orbits)[source]¶ Set and validate orbital parameters contain all required values.
Sets self.orbits and self.orb_format. If objid is not present in orbits, a sequential series of integers will be used. If H is not present in orbits, a default value of 20 will be used. If g is not present in orbits, a default value of 0.15 will be used. If sed_filename is not present in orbits, either C or S type will be assigned, according to the semi-major axis value.
- Parameters
orbits (pandas.DataFrame, pandas.Series or numpy.ndarray) – Array-like object containing orbital parameter information.
-
updateOrbits
(neworb)[source]¶ Update existing orbits with new values, leaving OrbitIds, H, g, and sed_filenames in place.
Example use: transform orbital parameters (using PyOrbEphemerides) and then replace original values. Example use: propagate orbital parameters (using PyOrbEphemerides) and then replace original values.
- Parameters
neworb (pandas.DataFrame) –
-
lsst.sims.maf.slicers.uniSlicer module¶
lsst.sims.maf.slicers.userPointsSlicer module¶
-
class
lsst.sims.maf.slicers.userPointsSlicer.
UserPointsSlicer
(ra, dec, lonCol='fieldRA', latCol='fieldDec', latLonDeg=True, verbose=True, badval=- 666, leafsize=100, radius=1.75, useCamera=False, rotSkyPosColName='rotSkyPos', mjdColName='observationStartMJD', chipNames='all')[source]¶ Bases:
lsst.sims.maf.slicers.baseSpatialSlicer.BaseSpatialSlicer
A spatial slicer that evaluates pointings overlapping user-provided list of points.
- Parameters
ra (list or numpy.ndarray) – User-selected RA points, in degrees. Stored internally in radians.
dec (list or numpy.ndarray) – User-selected Dec points, in degrees. Stored internally in radians.
lonCol (str, optional) – Name of the longitude (RA equivalent) column to use from the input data. Default fieldRA
latCol (str, optional) – Name of the latitude (Dec equivalent) column to use from the input data. Default fieldDec
latLonDeg (bool, optional) – Flag indicating whether the lon and lat values will be in degrees (True) or radians (False). Default True (appropriate for opsim v4).
verbose (boolean, optional) – Flag to indicate whether or not to write additional information to stdout during runtime. Default True.
badval (float, optional) – Bad value flag, relevant for plotting. Default -666.
leafsize (int, optional) – Leafsize value for kdtree. Default 100.
radius (float, optional) – Radius for matching in the kdtree. Equivalent to the radius of the FOV. Degrees. Default 1.75.
useCamera (boolean, optional) – Flag to indicate whether to use the LSST camera footprint or not. Default False.
rotSkyPosColName (str, optional) – Name of the rotSkyPos column in the input data. Only used if useCamera is True. Describes the orientation of the camera orientation compared to the sky. Default rotSkyPos.
mjdColName (str, optional) – Name of the exposure time column. Only used if useCamera is True. Default observationStartMJD.
chipNames (array-like, optional) – List of chips to accept, if useCamera is True. This lets users turn ‘on’ only a subset of chips. Default ‘all’ - this uses all chips in the camera.