The Mockimgs
module¶
The Mockimgs
modules provides tools to generate mock observation images from N-body models.
Designing an instrument¶
In its current form, an instrument is composed of a telescope, a ccd and a filter:
from astropy import units as u
from pNbody.Mockimgs import telescope, ccd, instruments, filters, obs_object
# define a telescope
mytelescope = telescope.Telescope(name="MyTelescope",focal=1500*u.mm)
# define a ccd
myccd = ccd.CCD(name="MyCCD",shape=[4096,4096],size=[41*u.mm,41*u.mm])
# define a filter
myfilter = filters.Filter("BastI_GAIA_G")
# define an instrument
myinstrument = instrument.Instrument(
name = "MyInstrument",
telescope = mytelescope,
ccd = myccd,
filter_type = myfilter,
)
# get info
myinstrument.info()
The last command will provide the following output:
name : MyInstrument
fov : [1.565987147680034 deg,1.565987147680034 deg]
telescope name : MyTelescope
ccd name : MyCCD
ccd fov : [1.565987147680034 deg,1.565987147680034 deg]
ccd shape : [4096,4096]
ccd xmin xmax : -0.013665815885465636 rad 0.013665815885465636 rad
ccd ymin ymax : -0.013665815885465636 rad 0.013665815885465636 rad
ccd xmin xmax (kpc): 1.0 kpc 1.0 kpc
ccd ymin ymax (kpc): 1.0 kpc 1.0 kpc
filter name : G
filter module name : <module 'pNbody.Mockimgs.Mag_BastI' from '../python3.9/site-packages/pNbody/Mockimgs/Mag_BastI.py'>
Note that the focal is used here to convert physical distances into angles:
\[\theta = \arctan( s / f)\]
where for example, \(\theta\) is the field of view in radian, \(s\) is the ccd size in some length units and \(f\) is the focal expressed in the same units.
Alternatively, pre-defined instruments (see the list below) can be directly loaded via
myinstrument = instruments["arrakihs_vis_G"]
myinstrument.info()
which gives:
name : arrakihs_vis_G
fov : [1.4997222222222222 deg,1.4997222222222222 deg]
telescope name : iSIM-170
ccd name : arrakihs_vis
ccd fov : [1.4997222222222222 deg,1.4997222222222222 deg]
ccd shape : [4096,4096]
ccd xmin xmax : -2699.5 arcsec 2699.5 arcsec
ccd ymin ymax : -2699.5 arcsec 2699.5 arcsec
ccd xmin xmax (kpc): 1.0 kpc 1.0 kpc
ccd ymin ymax (kpc): 1.0 kpc 1.0 kpc
filter name : G
filter module name : <module 'pNbody.Mockimgs.Mag_BastI' from '/home/revaz/.local/lib/python3.9/site-packages/pNbody/Mockimgs/Mag_BastI.py'>
Attach and object¶
To generate an image the instrument must be attached to an object (e.g. galaxy):
# create the object to observe ("snapshot.hdf5" is an N-body model)
galaxy = obs_object.Object(filename="snapshot.hdf5",distance=35*u.Mpc)
# add the object to the instrument
myinstrument.add_object(galaxy)
Compute flux and surface brightness maps¶
If the N-body model snapshot.hdf5
contains all necessary informations (see Setup the N-body model),
flux and surface brightness maps can be computed:
# set object for exposure
myinstrument.set_object_for_exposure()
# compute maps
FluxMap = myinstrument.ComputeFluxMap()
SBMap = myinstrument.SurfaceBrightnessMap()
# save a fits image
myinstrument.saveFitsImage(SBMap,"image.fits")
SSP Magnitudes in different filters¶
The Mockimgs
modules provides stellar magnitudes for single stellar populations (SSP) in different filters.
Those are provided by the Mockimgs.filters
.
Below is an example showing how to use it:
from pNbody.Mockimgs import filters
# get the list of all pNbody filters
filters.List()
# set a filter from the previous list
Filter = filters.Filter("BastI_GAIA_G_RP")
# or from a file
Filter = filters.Filter("pNbody/config/opt/filters/BPASS230_EuclidNISPY.hdf")
Mass = 1e11 # in Msol, the initial SSP mass
Age = 13 # in Gyr
MH = -2 # in log10([M/H])
# compute the magnitude
mag = Filter.Magnitudes(Mass,Age,MH)
print(mag)
-21.824502416627993
A list of available filters is provided here: List of filters.
SSP Bolometric luminosities¶
Bolometric luminosities can be accessed via the luminosities
module
from pNbody.Mockimgs import luminosities
# set a "filter"
filter_name = "BastI_L"
LM = luminosities.LuminosityModel(filter_name)
Mass = 1e11 # in Msol, the initial SSP mass
Age = 13 # in Gyr
MH = -2 # in log10([M/H])
# compute the magnitude
L = LM.Luminosities(Mass,Age,MH)
print(L)
27352237694.01902
Note
In the case where the initial SSP mass is not available, it is possible to provide the current SSP mass
and inform pNbody
that it must apply a correction, using the current_mass
option
Mag = Filter.Magnitudes(Mass,Age,MH,current_mass=True)
L = LM.Luminosities(Mass,Age,MH,current_mass=True)
The correction is done using a table stored in the hdf5 files which gives, for each age and metallicity the value of the current SSP mass compared to the initial value.
Pre-defined instruments¶
List of available filters¶
See List of filters
pNbody.Mockimgs
classes¶
The instrument
class¶
- class pNbody.Mockimgs.instrument.Instrument(name=None, telescope=None, ccd=None, ifu=None, mos=None, filter_type=None, obs_object=None, mapping_method=None)¶
Define the Instrument class. An instrument is composed of a telescope, a detector (ccd) and a filter.
- ComputeFluxMap(magFromField=None)¶
Compute a luminosity of flux map. Here the nbody model is supposed to be in the right units
It return a 2D matrix containing either a fluxes or masses or luminosities
If magFromField is not none, the magnitude is not computed but taken form the nbody object with the field name magFromField.
- FilterFluxMap(filter_opts={'mode': None})¶
Filter the flux map
FluxMas : the flux map (numpy.array) filter_opts : filter options
- SurfaceBrightnessMap()¶
Convert a mass/luminosity/flux to a surface brightness map
- add_object(obj)¶
add an object to observe
- change_fov(fov)¶
Use only a subset of the available fov
fov : new field of view
- draw(ax, mode=None, unit=None)¶
draw using matplotlib
- fov()¶
get the instrument (ccd) field of view
- getFluxMap()¶
return the current flux map
- getFocal()¶
return the telescope focal if defined
- getObject()¶
return the object
- getSurfaceBrightnessMap()¶
return the current surface brightness map
- get_arcsec_to_pixel_factors()¶
return the factor to convert arcsec to pixels
- get_parameters(fmt=None)¶
return a dictionary containing usefull parameters
- info()¶
get info on the instrument
- saveFitsImage(image, filename, units=None, comments=None, compress=True)¶
save a fits image
image : name of the image filename : output file name units : pixel units (str) comments : additional comments (list)
- set_filter(filtername)¶
set the filter
- filternamename of the filter, e.g.“BastI_GAIA_G”
of filter object , e.g. : filters.Filter(“BastI_GAIA_G”)
- set_mapping_method(mapping_method)¶
Set the kernel used to map particles to the CCD.
mapping_method: ‘gauss_old’, ‘gauss’, ‘spline’
- set_object_for_exposure(getAgeMH=False)¶
Open the N-body model and eventually compute missing quantities. Rotate according to a given los. Scale according to a given distance.
The telescope
class¶
- class pNbody.Mockimgs.telescope.Telescope(name=None, focal=None)¶
Define the telescope class.
- get_focal()¶
return the telescope focal
- get_parameters(fmt=None)¶
return a dictionary containing usefull parameters
- info()¶
give info on the telescope
- set_focal(focal)¶
set the telescope focal
focal : an astropy quantity in unit length
The ccd
class¶
- class pNbody.Mockimgs.ccd.CCD(name, shape=None, size=None, pixel_size=None, fov=None, focal=None)¶
Define the CCD class.
- MapUsingGauss(pos, amp, rsp)¶
create a map (using a gaussian kernel) from points with coordinates = pos : in arsec amplitude = amp size = rsp
- MapUsingSpline(pos, amp, rsp)¶
create a map (using a gaussian kernel) from points with coordinates = pos : in arsec amplitude = amp size = rsp
- change_fov(fov)¶
change the field of view, by modifying the shape
- defineSensibleZones()¶
set the sensible zone in size units
self.sZone is a list of tupple
- draw(ax, mode=None, unit=None)¶
draw the detector
- get_arcsec_to_pixel_factors()¶
return the factor to convert arcsec to pixels
- get_parameters(fmt=None)¶
return a dictionary containing usefull parameters
- info()¶
gives info on the ccd
- init()¶
Initialize all necessary quantities
- nx()¶
return the number of pixels along the x axis
- ny()¶
return the number of pixels along the y axis
- set_extension()¶
set the ccd extension in angle
- set_extension_in_spatial_length(object_distance)¶
set the ccd extension in spatial length
distance : distance to the observed object
- set_focal(focal)¶
set the focal and ignore the the fov
- set_focal_from_fov()¶
set the focal from the field of view and size
- set_fov_from_focal()¶
set the field of view from the focal
- set_pixel_area()¶
set the pixel area needs pixel_fov to be defined
- set_pixel_fov()¶
set the pixel field of view needs the focal to be defined
- set_pixel_size()¶
set size from pixel_size and shape
- set_size()¶
set size from pixel_size and shape
- set_sizeOrpixel_size()¶
set pixel_size from size or set size from pixel_size
The filters
module¶
- pNbody.Mockimgs.filters.Filter(name, filter_type=None)¶
the filter class builder
- class pNbody.Mockimgs.filters.FilterCl(filter_filename)¶
Define a filter class that describe a filter from its transmission curve.
filters are characterised by their transmission curve given as a function of the wavelength (by default in angstrom)
- SEDmultiply(SED, wavelength)¶
multipy a given SED with the current filter
SED : the SED in arbitrary units wavelength : wavelength in angstrom
- get_response()¶
get the filter response
- get_wavelength(units=Unit('Angstrom'))¶
get the wavelength
- init()¶
init a set of useful quantities defining the properties of the filter
- read()¶
Read the file that contains the filter information and return the filter transmission as a function of the wavelength in Angstrom
- class pNbody.Mockimgs.filters.FilterClass(filename, filter_name=None)¶
the filter class
- Read()¶
read the file and create a table
- get_parameters(fmt=None)¶
return a dictionary containing useful parameters
- class pNbody.Mockimgs.filters.FilterGAEA(filename, filter_name=None)¶
the filter GAEA class
- Magnitudes(mass, age, Z)¶
Ages are provided in Gyrs Mass in Msol
- Read()¶
special read for FilterGAEA
- class pNbody.Mockimgs.filters.FilterHDF5(filename, filter_name=None)¶
the Filter HDF5 class
- Magnitudes(mass, age, MH, current_mass=False)¶
mass : initial SSP mass in Msol age : SSP age in Gyr MH : SSP metallicity in log10(M/M0)
- current_massif True, assume the mass to be the final SSP mass
and convert it to the initial SSP mass default=False
: if False, assume the mass is the initial SSP mass
- Read()¶
special read for FilterPKL
- SSPMassRatio(age, MH)¶
Return the ratio of the initial SSP to the final SSP mass
- class pNbody.Mockimgs.filters.FilterPKL(filename, filter_name=None)¶
the Filter PKL class
- Magnitudes(mass, age, Z)¶
Ages are provided in Gyrs
- Read()¶
special read for FilterPKL
- pNbody.Mockimgs.filters.List()¶
return a list of current available filters
- pNbody.Mockimgs.filters.list()¶
return a list of current available filters
The obs_object
class¶
- class pNbody.Mockimgs.obs_object.Object(name=None, filename=None, unitsfile=None, nb=None, adjust_coords='remove_identical', distance=<Quantity 1. Mpc>, unit=Unit("kpc"), ref='universe', los=<pNbody.Mockimgs.los.LineOfSights object>, rsp_opts={'fac': 1, 'mode': None}, ftype='gh5', verbose=0, local=False, status='old')¶
Define an object “object”. the latter derives from the Nbody class
- ScaleSmoothingLength(x, opt)¶
Scale the smoothing length. opt : a dictionary
- align(pos, axis, axis_ref=[0, 0, 1])¶
rotate the model with a rotation that align axis with axis_ref
pos : a set of positions axis : the axis axis_ref : the reference axis
- draw(ax, mode=None, unit=None)¶
draw in matplotlib do a simple scatter plot
- getDistance()¶
return the object distance
- getFocal()¶
if defined, get the instrument focal
- get_parameters(fmt=None)¶
return a dictionary containing usefull parameters
- info()¶
give info on the telescope
- open(getAgeMH=False)¶
Open the N-body model and eventually compute missing quantities.
- set_for_exposure()¶
Set the object for exposure. Rotate according to a given los. Scale according to a given distance. Initial values of the model are modified.
- set_los(los)¶
set line of sights
- set_rsp_opts(rsp_opts)¶
set rsp options
- x(mode=None, unit=None)¶
get x coordinates of the object we assume that x is in kpc
- y(mode=None, unit=None)¶
get y coordinates of the object we assume that y is in kpc
- class pNbody.Mockimgs.obs_object.oldObject(name=None, filename=None, unitsfile=None, nb=None, adjust_coords='remove_identical', distance=<Quantity 1. Mpc>, unit=Unit("kpc"), ref='universe', los=<pNbody.Mockimgs.los.LineOfSights object>, rsp_opts={'fac': 1, 'mode': None})¶
Define the object class.
- ScaleSmoothingLength(x, opt)¶
Scale the smoothing length. opt : a dictionary
- align(pos, axis, axis_ref=[0, 0, 1])¶
rotate the model with a rotation that align axis with axis_ref
pos : a set of positions axis : the axis axis_ref : the reference axis
- draw(ax, mode=None, unit=None)¶
draw in matplotlib do a simple scatter plot
- getDistance()¶
return the object distance
- getFocal()¶
if defined, get the instrument focal
- get_parameters(fmt=None)¶
return a dictionary containing usefull parameters
- info()¶
give info on the telescope
- open()¶
Open the N-body model and eventually compute missing quantities.
- set_for_exposure()¶
Set the object for exposure. Rotate according to a given los. Scale according to a given distance. Initial values of the model are modified.
- set_los(los)¶
set line of sights
- set_rsp_opts(rsp_opts)¶
set rsp options
- x(mode=None, unit=None)¶
get x coordinates of the object we assume that x is in kpc
- y(mode=None, unit=None)¶
get y coordinates of the object we assume that y is in kpc