Spectral classes and functions
Classes and functions for managing and processing spectra
Core Functions
Core functions for spectral data.
- dysh.spectra.core.average(data, axis=0, weights=None)[source]
Average a group of spectra or scans.
- Parameters:
- data
ndarray
The spectral data, typically with shape (nspect,nchan).
- axisint
The axis over which to average the data. Default axis=0 will return the average spectrum if shape is (nspect,nchan)
- weights
ndarray
The weights to use in averaging. These might typically be system temperature based. The weights array must be the length of the axis over which the average is taken. Default: None will use equal weights
- data
- Returns:
- average
ndarray
The average along the input axis
- average
- dysh.spectra.core.baseline(spectrum, order, exclude=None, **kwargs)[source]
Fit a baseline for a spectrum
- Parameters:
- spectrum
Spectrum
The input spectrum
- orderint
The order of the polynomial series, a.k.a. baseline order
- excludelist of 2-tuples of int or
Quantity
, orSpectralRegion
List of region(s) to exclude from the fit. The tuple(s) represent a range in the form [lower,upper], inclusive. Examples:
One channel-based region:
>>> [11,51]
Two channel-based regions:
>>> [(11,51),(99,123)]
One
Quantity
region:>>> [110.198*u.GHz,110.204*u.GHz].
One compound
SpectralRegion
:>>> SpectralRegion([(110.198*u.GHz,110.204*u.GHz),(110.196*u.GHz,110.197*u.GHz)]). Default: no exclude region
- modelstr
One of ‘polynomial’ or ‘chebyshev’, Default: ‘polynomial’
- fitter
_FitterMeta
The fitter to use. Default:
LinearLSQFitter
(withcalc_uncertaintes=True
). Be care when choosing a different fitter to be sure it is optimized for this problem.
- spectrum
- Returns:
- modelslist of
Model
The list of models that contain the fitted model parameters. See
fit_continuum
.
- modelslist of
- dysh.spectra.core.exclude_to_region(exclude, refspec, fix_exclude=False)[source]
Convert an exclude list to a list of ~specutuls.SpectralRegion.
- Parameters:
- excludelist of 2-tuples of int or
Quantity
, orSpectralRegion
List of region(s) to exclude from the fit. The tuple(s) represent a range in the form [lower,upper], inclusive. Examples:
One channel-based region:
>>> [11,51]
Two channel-based regions:
>>> [(11,51),(99,123)]
One
Quantity
region:>>> [110.198*u.GHz,110.204*u.GHz].
One compound
SpectralRegion
:>>> SpectralRegion([(110.198*u.GHz,110.204*u.GHz),(110.196*u.GHz,110.197*u.GHz)]).
- refspec: `~spectra.spectrum.Spectrum`
The reference spectrum whose spectral axis will be used when converting between exclude and axis units (e.g. channels to GHz).
- fix_exclude: bool
If True, fix exclude regions that are out of bounds of the specctral axis to be within the spectral axis. Default:False
- excludelist of 2-tuples of int or
- Returns:
- regionlistlist of
SpectralRegion
A list of
SpectralRegion
corresponding toexclude
with units of therefspec.spectral_axis
.
- regionlistlist of
- dysh.spectra.core.mean_tsys(calon, caloff, tcal, mode=0, fedge=0.1, nedge=None)[source]
Get the system temperature from the neighboring calon and caloff, which reflect the state of the noise diode. We define an extra way to set the edge size, nedge, if you prefer to use number of edge channels instead of the inverse fraction. This implementation recreates GBTIDL’s
dcmeantsys
.- Parameters:
- calon
ndarray
-like ON calibration
- caloff
ndarray
-like OFF calibration
- tcal
ndarray
-like calibration temperature
- modeint
mode=0 Do the mean before the division mode=1 Do the mean after the division
- fedgefloat
Fraction of edge channels to exclude at each end, a number between 0 and 1. Default: 0.1, meaning the central 80% bandwidth is used
- nedgeint
Number of edge channels to exclude. Default: None, meaning use
fedge
- calon
- Returns:
- meanTsys
ndarray
-like The mean system temperature
- meanTsys
- dysh.spectra.core.region_to_axis_indices(region, refspec)[source]
- Parameters:
- region
SpectralRegion
- refspec: `~spectra.spectrum.Spectrum`
The reference spectrum whose spectral axis will be used when converting between exclude and axis units (e.g., channels to GHz).
- region
- Returns:
- indices2-tuple of int
The array indices in
refspec
corresponding toregion.bounds
- dysh.spectra.core.sq_weighted_avg(a, axis=0, weights=None)[source]
Compute the mean square weighted average of an array (2nd moment).
\(v = \sqrt{\frac{\sum_i{w_i~a_i^{2}}}{\sum_i{w_i}}}\)
- Parameters:
- Returns:
- average
ndarray
The square root of the squared average of a along the input axis.
- average
- dysh.spectra.core.tsys_weight(exposure, delta_freq, tsys)[source]
- Compute the system temperature based weight(s) using the expression:
\(w = t_{exp} \times \delta_\nu / T_{sys}^2,\)
If
exposure
,delta_freq
, ortsys
parameters are given asQuantity
, they will be converted to seconds, Hz, K, respectively for the calculation. (Not that this really matters since weights are relative to each other)>>> import astropy.units as u >>> [3*u.s, 4*u.s] ### WRONG
Rather, if using
Quantity
, they have to beQuantity
objects, e.g.,>>> [3, 4]*u.s ### RIGHT!
Spectrum
The Spectrum class to contain and manipulate spectra.
- class dysh.spectra.spectrum.Spectrum(*args, **kwargs)[source]
Bases:
Spectrum1D
This class contains a spectrum and its attributes. It is built on
Spectrum1D
with added attributes like baseline model. Note thatSpectrum1D
can contain multiple spectra but we probably will not use that because the restriction that it can have only one spectral axis conflicts with slight Doppler shifts. SeeSpectrum1D
for the instantiation arguments.- Attributes:
array_axis_physical_types
Returns the WCS physical types that vary along each array axis.
baseline_model
Returns the computed baseline model or None if it has not yet been computed.
- bin_edges
combined_wcs
The WCS transform for the NDCube, including the coordinates specified in
.extra_coords
.data
ndarray
-like : The stored dataset.- dimensions
doppler_convention
String representation of the velocity (Doppler) convention
energy
The energy of the spectral axis as a
Quantity
in units of eV.equivalencies
Get the spectral axis equivalencies that can be used in converting the axis
- exclude_regions
extra_coords
Coordinates not described by
NDCubeABC.wcs
which vary along one or more axes.flux
Converts the stored data and unit information into a quantity.
frequency
The
spectral_axis
as aQuantity
in units of GHzglobal_coords
Coordinate metadata which applies to the whole cube.
mask
any type : Mask for the dataset, if any.
- meta
dict
-like : Additional meta information about the dataset.observer
Returns ——- observer :
BaseCoordinateFrame
or derivative The coordinate frame of the observer if present.- obstime
photon_flux
The flux density of photons as a
Quantity
, in units of- plotter
psf
Image representation of the PSF for the dataset.
radial_velocity
The radial velocity(s) of the objects represented by this spectrum.
redshift
The redshift(s) of the objects represented by this spectrum.
- rest_value
- shape
spectral_axis
Returns the SpectralCoord object.
- spectral_axis_direction
spectral_axis_unit
Deprecated since version v1.1.
spectral_wcs
Returns the spectral axes of the WCS
- target
uncertainty
any type : Uncertainty in the dataset, if any.
unit
Unit
: Unit for the dataset, if any.velocity
Converts the spectral axis array to the given velocity space unit given the rest value.
velocity_convention
Returns the velocity convention
velocity_frame
String representation of the velocity frame
wavelength
The
spectral_axis
as aQuantity
in units of Angstromswcs
any type : A world coordinate system (WCS) for the dataset, if any.
Methods
add
(operand[, operand2])Performs addition by evaluating
self
+operand
.axis_velocity
([unit])Get the spectral axis in velocity units.
axis_world_coords
(*axes[, pixel_corners, wcs])Returns objects representing the world coordinates of pixel centers for a desired axes.
axis_world_coords_values
(*axes[, ...])Returns the world coordinate values of all pixels for desired axes.
baseline
(degree[, exclude])Compute and optionally remove a baseline.
bshow
()Show the baseline model
collapse
(method[, axis])Collapse the flux array given a method.
crop
(*points[, wcs])Crop using real world coordinates.
crop_by_values
(*points[, units, wcs])Crop using real world coordinates.
divide
(operand[, operand2])Performs division by evaluating
self
/operand
.explode_along_axis
(axis)Separates slices of NDCubes along a given axis into an NDCubeSequence of (N-1)DCubes.
make_spectrum
(data, meta[, use_wcs, ...])Factory method to create a Spectrum object from a data and header.
multiply
(operand[, operand2])Performs multiplication by evaluating
self
*operand
.new_flux_unit
(unit[, equivalencies, ...])Converts the flux data to the specified unit.
plot
(**kwargs)Visualize the
NDCube
.rebin
(bin_shape[, operation, ...])Downsample array by combining contiguous pixels into bins.
reproject_to
(target_wcs[, algorithm, ...])Reprojects the instance to the coordinates described by another WCS object.
savefig
(file, **kwargs)Save the plot
set_convention
(doppler_convention)Set the velocity convention of this Spectrum.
set_frame
(toframe)Set the sky coordinate and doppler tracking reference frame of this Spectrum.
set_radial_velocity_to
(radial_velocity)This sets the radial velocity of the spectrum to be
radial_velocity
without changing the values of thespectral_axis
.set_redshift_to
(redshift)This sets the redshift of the spectrum to be
redshift
without changing the values of thespectral_axis
.shift_spectrum_to
(*[, redshift, radial_velocity])This shifts in-place the values of the
spectral_axis
, given either a redshift or radial velocity.stats
()Compute some statistics of this
Spectrum
.subtract
(operand[, operand2])Performs subtraction by evaluating
self
-operand
.to
(new_unit, **kwargs)Convert instance to another unit.
velocity_axis_to
([unit, toframe, ...])- Parameters:
with_frame
(toframe)Return a copy of this Spectrum with a new coordinate reference frame.
with_spectral_unit
(unit[, ...])Returns a new spectrum with a different spectral axis unit.
with_velocity_convention
(doppler_convention)Returns a copy of this Spectrum with the input velocity convention.
get_velocity_shift_to
list_to_spectral_region
max
mean
median
min
read
smooth
sum
write
- axis_velocity(unit=Unit('km / s'))[source]
Get the spectral axis in velocity units. Note: This is not the same as
Spectrum.velocity
, which includes the source radial velocity.
- baseline(degree, exclude=None, **kwargs)[source]
Compute and optionally remove a baseline. The model for the baseline can be either a 1D polynomial model or a 1D Chebyshev polynomial of the first kind. The code uses
astropy.modeling
andastropy.fitter
to compute the baseline. See the documentation for those modules. This method will set thebaseline_model
attribute to the fitted model function which can be evaluated over a domain.- Parameters:
- degreeint
The degree of the polynomial series, a.k.a. baseline order
- excludelist of 2-tuples of int or ~astropy.units.Quantity, or ~specutils.SpectralRegion
List of region(s) to exclude from the fit. The tuple(s) represent a range in the form [lower,upper], inclusive. In channel units.
Examples:
One channel-based region: [11,51]
Two channel-based regions: [(11,51),(99,123)].
One ~astropy.units.Quantity region: [110.198*u.GHz,110.204*u.GHz].
One compound
SpectralRegion
: SpectralRegion([(110.198*u.GHz,110.204*u.GHz),(110.196*u.GHz,110.197*u.GHz)]).Default: no exclude region
- modelstr
One of ‘polynomial’ or ‘chebyshev’, Default: ‘polynomial’
- fitter
_FitterMeta
The fitter to use. Default:
LinearLSQFitter
(withcalc_uncertaintes=True
). Be care when choosing a different fitter to be sure it is optimized for this problem.- removebool
If True, the baseline is removed from the spectrum. Default: False
- property baseline_model
Returns the computed baseline model or None if it has not yet been computed.
- property doppler_convention
String representation of the velocity (Doppler) convention
- property equivalencies
Get the spectral axis equivalencies that can be used in converting the axis between km/s and frequency or wavelength
- property exclude_regions
- classmethod make_spectrum(data, meta, use_wcs=True, observer_location=None)[source]
Factory method to create a Spectrum object from a data and header.
- Parameters:
- data
ndarray
The data array. See
Spectrum1D
- metadict
The metadata, typically derived from an SDFITS header. Required items in
meta
are ‘CTYPE[123]’,’CRVAL[123]’, ‘CUNIT[123]’, ‘VELOCITY’, ‘EQUINOX’, ‘RADESYS’- use_wcsbool
If True, create a WCS object from
meta
- observer_location
EarthLocation
Location of the observatory. See
Observatory
. This will be transformed toITRS
using the time of observation DATE-OBS or MJD-OBS inmeta
.
- data
- Returns:
- spectrum
Spectrum
The spectrum object
- spectrum
- property observer
- Returns:
- observer
BaseCoordinateFrame
or derivative - The coordinate frame of the observer if present.
- observer
- property obstime
- plot(**kwargs)[source]
Visualize the
NDCube
.- Parameters:
- axes: `~astropy.visualization.wcsaxes.WCSAxes` or None:, optional
The axes to plot onto. If None the current axes will be used.
- plot_axes: `list`, optional
A list of length equal to the number of pixel dimensions in array axis order. This list selects which cube axes are displayed on which plot axes. For an image plot this list should contain
'x'
and'y'
for the plot axes andNone
for all the other elements. For a line plot it should only contain'x'
andNone
for all the other elements.- axes_unit: `list`, optional
A list of length equal to the number of world dimensions specifying the units of each axis, or
None
to use the default unit for that axis.- axes_coordinates: `list`, optional
A list of length equal to the number of pixel dimensions. For each axis the value of the list should either be a string giving the world axis type or
None
to use the default axis from the WCS.- data_unit: `astropy.units.Unit`
The data is changed to the unit given or the
NDCube.unit
if not given.- wcs: `astropy.wcs.wcsapi.BaseHighLevelWCS`
The WCS object to define the coordinates of the plot axes.
- kwargs
Additional keyword arguments are given to the underlying plotting infrastructure which depends on the dimensionality of the data and whether 1 or 2 plot_axes are defined: - Animations:
mpl_animators.ArrayAnimatorWCS
- Static 2-D images:matplotlib.pyplot.imshow
- Static 1-D line plots:matplotlib.pyplot.plot
- property plotter
- set_convention(doppler_convention)[source]
Set the velocity convention of this Spectrum. The spectral axis of this Spectrum will be replaced with a new spectral axis with the input velocity convention. The header ‘VELDEF’ value will be changed accordingly.
To make a copy of this Spectrum with a new velocity convention instead, use
with_velocity_convention
.- Parameters:
- doppler_convention - str
The velocity convention, one of ‘radio’, ‘optical’, ‘relativistic’
- set_frame(toframe)[source]
Set the sky coordinate and doppler tracking reference frame of this Spectrum. The header ‘CTYPE1’ will be changed accordingly.
To make a copy of this Spectrum with new coordinate referece frmae instead, use
with_frame
.- Parameters:
- toframe - str
The coordinate reference frame identifying string, as used by astropy, e.g. ‘hcrs’, ‘icrs’, etc.
- stats()[source]
Compute some statistics of this
Spectrum
. The mean, rms, data minimum and data maximum are calculated. Note this works with slicing, so, e.g.,myspectrum[45:153].stats()
will return the statistics of the slice.- Returns:
- statstuple
Tuple consisting of (mean,rms,datamin,datamax)
- property target
- velocity_axis_to(unit=Unit('km / s'), toframe=None, doppler_convention=None)[source]
- Parameters:
- unit
Quantity
or str that can be converted to Quantity The unit to which the axis is to be converted
- toframestr
The coordinate frame to convert to, e.g. ‘hcrs’, ‘icrs’
- doppler_conventionstr
The Doppler velocity covention to use, one of ‘optical’, ‘radio’, or ‘rest’
- unit
- Returns:
- velocity
Quantity
The converted spectral axis velocity
- velocity
- property velocity_frame
String representation of the velocity frame
- with_frame(toframe)[source]
Return a copy of this Spectrum with a new coordinate reference frame.
- Parameters:
- toframe - str
The coordinate reference frame identifying string, as used by astropy, e.g. ‘hcrs’, ‘icrs’, etc.
- Returns:
- spectrum
Spectrum
A new Spectrum object
- spectrum
- with_velocity_convention(doppler_convention)[source]
Returns a copy of this Spectrum with the input velocity convention. The header ‘VELDEF’ value will be changed accordingly.
- Parameters:
- doppler_convention - str
The velocity convention, one of ‘radio’, ‘optical’, ‘relativistic’
- Returns:
- spectrum
Spectrum
A new Spectrum object
- spectrum
Scan
The classes that define various types of Scan and their calibration methods.
- class dysh.spectra.scan.FSScan(gbtfits, scan, sigrows, calrows, bintable, calibrate=True, fold=True, use_sig=True, observer_location=<EarthLocation (882593.9465029, -4924896.36541728, 3943748.74743984) m>, debug=False)[source]
Bases:
ScanMixin
GBT specific version of Frequency Switch Scan
- Parameters:
- gbtfits
GBFITSLoad
input GBFITSLoad object
- scanint
scan number that contains integrations with a series of sig/ref and calon/caloff states
- sigrows: dict
dictionary containing with keys ‘ON’ and ‘OFF’ containing list of rows in
sdfits
corresponding to sig=T (ON) and sig=F (OFF) integrations.- calrowsdict
dictionary containing with keys ‘ON’ and ‘OFF’ containing list of rows in
sdfits
corresponding to cal=T (ON) and cal=F (OFF) integrations.- bintableint
the index for BINTABLE in
sdfits
containing the scans- calibrate: bool
whether or not to calibrate the data. If true, data will be calibrated as TSYS*(ON-OFF)/OFF. Default: True
- fold: bool
whether or not to fold the spectrum. Default: True
- use_sigbool
whether to use the sig as the sig, or the ref as the sig. Default: True
- observer_location
EarthLocation
Location of the observatory. See
Observatory
. This will be transformed toITRS
using the time of observation DATE-OBS or MJD-OBS in the SDFITS header. The default is the location of the GBT.
- gbtfits
- Attributes:
delta_freq
Get the array of channel frequency width
exposure
Get the array of exposure (integration) times for FSscan
fdnum
The feed number
folded
Has the FSscan been folded?
ifnum
The IF number
nchan
The number of channels in this scan
npol
The number of polarizations in this Scan
nrows
The number of rows in this Scan
pols
The polarization number(s)
scan
The scan number
tsys
The system temperature array.
Methods
calibrate
(**kwargs)Frequency switch calibration, following equations .
calibrated
(i)Return the calibrated Spectrum of this FSscan
finalspectrum
([weights])Average all times and polarizations in this Scan
polaverage
([weights])Average all polarizations in this Scan
timeaverage
([weights])Compute the time-averaged spectrum for this set of FSscans.
- calibrate(**kwargs)[source]
Frequency switch calibration, following equations … fold=True or fold=False is required
- calibrated(i)[source]
Return the calibrated Spectrum of this FSscan
- Parameters:
- iint
The index into the calibrated array
- Returns:
- spectrum
Spectrum
- spectrum
- property delta_freq
Get the array of channel frequency width
df = [ 0.5*(df_ref_on + df_ref_off) + 0.5*(df_sig_on + df_sig_off) ] / 2
- Returns:
- delta_freq: ~numpy.ndarray
The channel frequency width in units of the CDELT1 keyword in the SDFITS header
- property exposure
Get the array of exposure (integration) times for FSscan
exposure = [ 0.5*(exp_ref_on + exp_ref_off) + 0.5*(exp_sig_on + exp_sig_off) ] / 2
- Returns:
- exposure~numpy.ndarray
The exposure time in units of the EXPOSURE keyword in the SDFITS header
- property folded
Has the FSscan been folded?
- Returns:
- boolean
True if the signal and reference integrations have been folded. False if not.
- class dysh.spectra.scan.PSScan(gbtfits, scans, scanrows, calrows, bintable, calibrate=True, observer_location=<EarthLocation (882593.9465029, -4924896.36541728, 3943748.74743984) m>)[source]
Bases:
ScanMixin
GBT specific version of Position Switch Scan. A position switch scan object has one IF, one feed, and one or more polarizations.
- Parameters:
- gbtfits
GBFITSLoad
input GBFITSLoad object
- scansdict
dictionary with keys ‘ON’ and ‘OFF’ containing unique list of ON (signal) and OFF (reference) scan numbers NOTE: there should be one ON and one OFF, a pair @todo ‘scans’ should become ‘scan’
- scanrowsdict
dictionary with keys ‘ON’ and ‘OFF’ containing the list of rows in
sdfits
corresponding to ON (signal) and OFF (reference) integrations- calrowsdict
dictionary containing with keys ‘ON’ and ‘OFF’ containing list of rows in
sdfits
corresponding to cal=T (ON) and cal=F (OFF) integrations.- bintableint
the index for BINTABLE in
sdfits
containing the scans- calibrate: bool
whether or not to calibrate the data. If true, data will be calibrated as TSYS*(ON-OFF)/OFF. Default: True
- observer_location
EarthLocation
Location of the observatory. See
Observatory
. This will be transformed toITRS
using the time of observation DATE-OBS or MJD-OBS in the SDFITS header. The default is the location of the GBT.
- gbtfits
- Attributes:
delta_freq
Get the array of channel frequency width
exposure
Get the array of exposure (integration) times
fdnum
The feed number
ifnum
The IF number
nchan
The number of channels in this scan
npol
The number of polarizations in this Scan
nrows
The number of rows in this Scan
pols
The polarization number(s)
scan
The scan number
scans
The dictionary of the ON and OFF scan numbers in the PSScan.
tsys
The system temperature array.
Methods
calibrate
(**kwargs)Position switch calibration, following equations 1 and 2 in the GBTIDL calibration manual
calibrated
(i)Return the calibrated Spectrum.
finalspectrum
([weights])Average all times and polarizations in this Scan
polaverage
([weights])Average all polarizations in this Scan
timeaverage
([weights])Compute the time-averaged spectrum for this set of scans.
- calibrate(**kwargs)[source]
Position switch calibration, following equations 1 and 2 in the GBTIDL calibration manual
- calibrated(i)[source]
Return the calibrated Spectrum.
- Parameters:
- iint
The index into the calibrated array
- Returns:
- spectrum
Spectrum
- spectrum
- property delta_freq
Get the array of channel frequency width
df = [ 0.5*(df_ref_on + df_ref_off) + 0.5*(df_sig_on + df_sig_off) ] / 2
- Returns:
- delta_freq: ~numpy.ndarray
The channel frequency width in units of the CDELT1 keyword in the SDFITS header
- property exposure
Get the array of exposure (integration) times
exposure = [ 0.5*(exp_ref_on + exp_ref_off) + 0.5*(exp_sig_on + exp_sig_off) ] / 2
- Returns:
- exposure~numpy.ndarray
The exposure time in units of the EXPOSURE keyword in the SDFITS header
- property scans
The dictionary of the ON and OFF scan numbers in the PSScan.
Returns scans : dict
The scan number dictionary
- class dysh.spectra.scan.ScanBlock(*args)[source]
-
- Attributes:
fdnum
The feed number
ifnum
The IF number
nchan
The number of channels in this scan
npol
The number of polarizations in this Scan
nrows
The number of rows in this Scan
pols
The polarization number(s)
scan
The scan number
Methods
append
(item)S.append(value) -- append value to the end of the sequence
calibrate
(**kwargs)Calibrate all scans in this ScanBlock
clear
()count
(value)extend
(other)S.extend(iterable) -- extend sequence by appending elements from the iterable
finalspectrum
([weights])Average all times and polarizations in all scans this ScanBlock
index
(value, [start, [stop]])Raises ValueError if the value is not present.
insert
(i, item)S.insert(index, value) -- insert value before index
polaverage
([weights])Average all polarizations in all scans in this ScanBlock
pop
([index])Raise IndexError if list is empty or index is out of range.
remove
(item)S.remove(value) -- remove first occurrence of value.
reverse
()S.reverse() -- reverse IN PLACE
timeaverage
([weights, mode])Compute the time-averaged spectrum for all scans in this ScanBlock.
copy
sort
- finalspectrum(weights='tsys')[source]
Average all times and polarizations in all scans this ScanBlock
- Parameters:
- weights: str
‘tsys’ or None. If ‘tsys’ the weight will be calculated as:
\(w = t_{exp} \times \delta\nu/T_{sys}^2\)
Default: ‘tsys’
- Returns
- ——-
- finalspectra: list of `~spectra.spectrum.Spectrum`
List of all the time- and polarization-averaged spectra
- polaverage(weights='tsys')[source]
Average all polarizations in all scans in this ScanBlock
- Parameters:
- weights: str
‘tsys’ or None. If ‘tsys’ the weight will be calculated as:
\(w = t_{exp} \times \delta\nu/T_{sys}^2\)
Default: ‘tsys’
- Returns
- ——-
- polaverage: list of `~spectra.spectrum.Spectrum`
List of all the polarization-averaged spectra
- timeaverage(weights='tsys', mode='old')[source]
Compute the time-averaged spectrum for all scans in this ScanBlock.
- Parameters:
- weights: str
‘tsys’ or None. If ‘tsys’ the weight will be calculated as:
\(w = t_{exp} \times \delta\nu/T_{sys}^2\)
Default: ‘tsys’
- Returns
- ——-
- timeaverage: list of `~spectra.spectrum.Spectrum`
List of all the time-averaged spectra
- class dysh.spectra.scan.ScanMixin[source]
Bases:
object
This class describes the common interface to all Scan classes. A Scan represents one IF, one feed, and one or more polarizations. Derived classes must implement
calibrate()
.- Attributes:
Methods
calibrate
(**kwargs)Calibrate the Scan data
finalspectrum
([weights])Average all times and polarizations in this Scan
polaverage
([weights])Average all polarizations in this Scan
timeaverage
([weights])Compute the time-averaged spectrum for this scan.
- property fdnum
The feed number
- Returns:
- int
The index of the Feed
- property ifnum
The IF number
- Returns:
- int
The index of the Intermediate Frequency
- property nchan
The number of channels in this scan
- Returns:
- int
The number of channels in this scan
- property npol
The number of polarizations in this Scan
- Returns:
- int
The number of polarizations in this Scan
- property nrows
The number of rows in this Scan
- Returns:
- int
The number of rows in this Scan
- property pols
The polarization number(s)
- Returns:
- list
The list of integer polarization number(s)
- property scan
The scan number
- Returns:
- int
The scan number of the integrations in the Scan object
- class dysh.spectra.scan.SubBeamNodScan(sigtp, reftp, fulltp=None, method='cycle', calibrate=True, observer_location=<EarthLocation (882593.9465029, -4924896.36541728, 3943748.74743984) m>, **kwargs)[source]
Bases:
ScanMixin
- Parameters:
- sigtp: list of ~spectra.scan.TPScan
Signal total power scans
- reftp: list ~spectra.scan.TPScan
Reference total power scans
- fulltp: ~spectra.scan.TPScan
A full (sig+ref) total power scans, used only for method=’scan’
- method: str
Method to use when processing. One of ‘cycle’ or ‘scan’. ‘cycle’ is more accurate and averages data in each SUBREF_STATE cycle. ‘scan’ reproduces GBTIDL’s snodka function which has been shown to be less accurate. Default:’cycle’
- calibrate: bool
Whether or not to calibrate the data.
- observer_location
EarthLocation
Location of the observatory. See
Observatory
. This will be transformed toITRS
using the time of observation DATE-OBS or MJD-OBS in the SDFITS header. The default is the location of the GBT.- weights: str
Weighting scheme to use when averaging the signal and reference scans ‘tsys’ or None. If ‘tsys’ the weight will be calculated as:
\(w = t_{exp} \times \delta\nu/T_{sys}^2\)
Default: ‘tsys’
- Attributes:
- delta_freq
- exposure
fdnum
The feed number
ifnum
The IF number
nchan
The number of channels in this scan
npol
The number of polarizations in this Scan
nrows
The number of rows in this Scan
pols
The polarization number(s)
scan
The scan number
- tsys
Methods
calibrate
(**kwargs)Calibrate the Scan data
finalspectrum
([weights])Average all times and polarizations in this Scan
polaverage
([weights])Average all polarizations in this Scan
timeaverage
([weights])Compute the time-averaged spectrum for this scan.
calibrated
- property delta_freq
- property exposure
- timeaverage(weights='tsys')[source]
Compute the time-averaged spectrum for this scan.
- Parameters:
- weights: str
‘tsys’ or None. If ‘tsys’ the weight will be calculated as:
\(w = t_{exp} \times \delta\nu/T_{sys}^2\)
Default: ‘tsys’
- Returns
- ——-
- spectrum
Spectrum
The time-averaged spectrum
- property tsys
- class dysh.spectra.scan.TPScan(gbtfits, scan, sigstate, calstate, scanrows, calrows, bintable, calibrate=True, observer_location=<EarthLocation (882593.9465029, -4924896.36541728, 3943748.74743984) m>)[source]
Bases:
ScanMixin
GBT specific version of Total Power Scan
- Parameters:
- gbtfits
GBFITSLoad
input GBFITSLoad object
- scan: int
scan number
- sigstatestr
one of ‘SIG’ or ‘REF’ to indicate if this is the signal or reference scan or ‘BOTH’ if it contains both
- calstatestr
one of ‘ON’ or ‘OFF’ to indicate the calibration state of this scan, or ‘BOTH’ if it contains both
- scanrowslist-like
the list of rows in
sdfits
corresponding to sig_state integrations- calrowsdict
dictionary containing with keys ‘ON’ and ‘OFF’ containing list of rows in
sdfits
corresponding to cal=T (ON) and cal=F (OFF) integrations forscan
- bintableint
the index for BINTABLE in
sdfits
containing the scans- calibrate: bool
whether or not to calibrate the data. If
True
, the data will be (calon - caloff)*0.5, otherwise it will be SDFITS row data. Default:True
- gbtfits
- Attributes:
- calstate
delta_freq
Get the array of channel frequency width
exposure
Get the array of exposure (integration) times
fdnum
The feed number
ifnum
The IF number
nchan
The number of channels in this scan
npol
The number of polarizations in this Scan
nrows
The number of rows in this Scan
pols
The polarization number(s)
scan
The scan number
- sigstate
tsys
The system temperature array.
Methods
calc_tsys
(**kwargs)Calculate the system temperature array
calibrate
(**kwargs)Calibrate the Scan data
finalspectrum
([weights])Average all times and polarizations in this Scan
polaverage
([weights])Average all polarizations in this Scan
timeaverage
([weights])Compute the time-averaged spectrum for this set of scans.
total_power
(i)Return the total power spectrum
tpmeta
- property calstate
- property delta_freq
Get the array of channel frequency width
df = 0.5*(df_ref_on + df_ref_off)
Note we only have access to the refon and refoff row indices so can’t use sig here. This is probably incorrect
- Returns:
- delta_freq:
ndarray
The channel frequency width in units of the CDELT1 keyword in the SDFITS header
- delta_freq:
- property exposure
Get the array of exposure (integration) times
exposure = 0.5*(exp_ref_on + exp_ref_off)
Note we only have access to the refon and refoff row indices so can’t use sig here. This is probably incorrect
- Returns:
- exposure
ndarray
The exposure time in units of the EXPOSURE keyword in the SDFITS header
- exposure
- property sigstate
- timeaverage(weights='tsys')[source]
Compute the time-averaged spectrum for this set of scans.
- Parameters:
- weights: str
‘tsys’ or None. If ‘tsys’ the weight will be calculated as:
\(w = t_{exp} \times \delta\nu/T_{sys}^2\)
Default: ‘tsys’
- Returns
- ——-
- spectrum
Spectrum
The time-averaged spectrum