:py:mod:`kinms_fitter` ====================== .. py:module:: kinms_fitter Submodules ---------- .. toctree:: :titlesonly: :maxdepth: 1 kinms_fitter/index.rst prior_funcs/index.rst sb_profs/index.rst transformClouds/index.rst velocity_profs/index.rst warp_funcs/index.rst Package Contents ---------------- Classes ~~~~~~~ .. autoapisummary:: kinms_fitter.kinms_fitter kinms_fitter.sb_profs kinms_fitter.velocity_profs kinms_fitter.warp_funcs Functions ~~~~~~~~~ .. autoapisummary:: kinms_fitter.transformClouds .. py:class:: kinms_fitter(filename, spatial_trim=None, spectral_trim=None, linefree_chans=[1, 5]) Wrapper for easily kinematically modelling datacubes (typically from ALMA or VLA) using KinMS. ... Changeable Attributes --------------------- xc_guess : float Guess for Kinematic/morphological centre of the disc. RA, in degrees. [Default = Centre of supplied cube after clipping] yc_guess : float Guess for Kinematic/morphological centre of the disc. Dec, in degrees. [Default = Centre of supplied cube after clipping] vsys_guess : float Guess for Systemic velocity of the gas, in the radio velocity convention. km/s. [Default = Velocity of central channel of supplied cube after clipping] totflux_guess : float Total flux guess. Units the same as those used in the cube (e.g. Jy km/s). [Default = total sum of the input cube] velDisp_guess : float Velocity dispersion guess. km/s. [Default = 8 km/s] pa_guess : float Optional position angle guess, if no warp specified. Degrees. [Default = 0 deg] inc_guess : float Optional Inclination angle guess, if no warp specified. Degrees. [Default = 45 deg] inc_profile : list of warp_profile objects Optional variation of the inclination with radius as described by warp_profile objects. [Default = None] pa_profile : list of warp_profile objects Optional variation of the inclination with radius as described by warp_profile objects. [Default = None] sb_profile : list of sb_profs objects Optional - Define the surface brightness profile with a list of sb_profs objects which will be evaluated in order. [Default = single exponential disc] expscale_guess : float Optional - defines the scalelength of the default single exponential disc used if sb_profile is not set. arcseconds. [Default = a fifth of the total extent of the cube] skySampClouds : ndarray Optional - Input cloud model made with the skySampler tool if you dont want to fit the gas distribution. vel_profile : list of vel_profs objects Define the circular velocity profile with a list of vel_profs objects. [Default = tilted rings] radial_motion : list of radial_motion objects Optional - Define the radial motion types with radial_motion objects (see base KinMS documentation). [Default: None] nrings : int Optional - Number of tilted rings to use, spaced a beamwidth apart. [Default = fill the total extent of the cube] vel_guess : float Optional - starting guess for the tilted ring asymptotic velocity. km/s [Default = 28percent of the velocity width coverved by the cube] xcent_range : two element ndarray Allowed range for fitting of the kinematic centre. RA in degrees [Default: minimum and maximum RA within the input cube] ycent_range : two element ndarray Allowed range for fitting of the kinematic centre. Dec in degrees [Default: minimum and maximum Dec within the input cube] vsys_range : two element ndarray Allowed range for fitting of the systemic velocity. km/s. [Default: minimum and maximum velocity contained within the input cube] totflux_range : two element ndarray Allowed range for fitting of the total flux. Units the same as input cube. [Default: Min=0 and Max=3x the total flux in the input cube] velDisp_range : two element ndarray Allowed range for fitting of the velocity dispersion. km/s. [Default: Min=0 and Max=50 km/s] pa_range : two element ndarray Optional - Allowed range for fitting of the position angle. Degrees. [Default: Min=0 and Max=360 degrees] inc_range : two element ndarray Optional - Allowed range for fitting of the inclination angle. Degrees. [Default: Min=0 and Max=89 degrees] expscale_range : two element ndarray Optional - Allowed range for fitting of the default exponential scale radius. arcseconds. [Default: Min=0 and Max=the size of your input cube] vel_range : two element ndarray Optional - Allowed range for fitting of the tilted ring velocities. km/s. [Default: Min=0 and Max=half the cube total velocity width] pa_prior : object Optional - supply an object defining the PA prior (see GAStimator documentation). [Default = boxcar between range limits] xc_prior : object Optional - supply an object defining the RA centre prior (see GAStimator documentation). [Default = boxcar between range limits] yc_prior : object Optional - supply an object defining the Dec centre prior (see GAStimator documentation). [Default = boxcar between range limits] vsys_prior : object Optional - supply an object defining the Vsys prior (see GAStimator documentation). [Default = boxcar between range limits] inc_prior : object Optional - supply an object defining the inc prior (see GAStimator documentation). [Default = boxcar between range limits] totflux_prior : object Optional - supply an object defining the totflux prior (see GAStimator documentation). [Default = boxcar between range limits] velDisp_prior : object Optional - supply an object defining the velocity dispersion prior (see GAStimator documentation). [Default = boxcar between range limits] initial_guesses : ndarray Override the initial guesses used for the fit. [Default = No override] Control Attributes ------------------ objname : string Name of the object you are fitting. [Default = read from FITS header OBJECT entry, or 'object'] nprocesses : int Number of processors to use for MCMC fitting. [Default = fall back on GAStimator default choices] niters : int Number of iterations to use when fitting. [Default=3000] output_initial_model : bool Save the initially generated model as a FITS file. [Default=False] pdf : bool Make PDF outputs. [Default = True] silent : bool Fit silently. [Default = False] show_corner : bool Show corner plots if MCMC fitting. [Default = True] nSamps : int Number of samples to use for KinMS model. [Default = 5e5] chi2_var_correct : bool Use chi-sqr variance correction. [Default = True] text_output : bool Save output text files. [Default = True] save_all_accepted : bool Save all the accepted guesses (allowing you to recreate corner plots etc). [Default = True] interactive : bool Use interactive mode for plots. [Default = True] show_plots : bool Show the diagnotic plots. [Default = True] output_cube_fileroot : string File root for output best fit datacube. [Default = False] lnlike_func : object Override the default likelihood function. [Default =. Dont override] tolerance : float Tolerance for simple fit. Smaller numbers are more stringent (longer runtime) [Default=0.1] Informational Attributes (dont change unless you know what you are doing) ------------------------------------------------------------------------- bunit : string Unit of datacube fluxes cellsize : float Cube cellsize in arcsec timetaken : float Time taken for the fit in seconds. chi2_correct_fac : float chi2 correction factor applied. mask_sum : 0 Number of detected pixels in the data mask. labels : ndarray of strings Label for each fitted quantity .. py:method:: colorbar(mappable, ticks=None) Add a colorbar to a given plot param mappable : matplotlib figure Figure which to add the colorbar too param ticks : ndarray Ticks to show on the colourbar. [Default = None] .. py:method:: read_in_a_cube(path) Reads in the datacube. .. py:method:: get_header_coord_arrays(hdr) Get coordinate arrays from a FITS header .. py:method:: rms_estimate(cube, chanstart, chanend) Estimate the RMS in the inner quarter of the datacube in a given channel range. .. py:method:: from_fits_history(hdr) Get beam parameters if they happen to be in history keywords. Stolen from radio_beam, with thanks! .. py:method:: read_primary_cube(cube) Wrapper method for reading in datacube. .. py:method:: setup_params() Setup the fit guesses, minima, maxima etc from the inputs. .. py:method:: clip_cube() Clip the input datacube to size. .. py:method:: model_simple(param, fileName='') Function to create a model from the given input parameters. .. py:method:: mcmc_fit(initial_guesses, labels, minimums, maximums, fixed, priors, precision) Function to run the MCMC fit. .. py:method:: simple_chi2(theargs, info) Likelihood function for the simple fit mode. .. py:method:: logo() Returns the KinMS logo. .. py:method:: simple_fit(initial_guesses, labels, minimums, maximums, fixed) Runs the simple fit. .. py:method:: plot(block=True, overcube=None, savepath=None, **kwargs) Makes the plots. .. py:method:: write_text(bestvals, errup, errdown, units, fixed, runtime, mode, errors_warnings='None', fname='KinMS_fitter_output.txt') Write output text files. .. py:method:: run(method='mcmc', justplot=False, savepath='./', **kwargs) Run the fit. :param method: Type of fit to run. :type method: string - choice of 'mcmc', 'simple' or 'both'. :param justplot: Just show the initial model plot without fitting. Helpful for setting initial guesses. :type justplot: bool :param savepath: Path to append to saved files. :type savepath: string :param **kwargs: any additional keywords to pass to other methods. .. py:class:: sb_profs Creates flexible surface brightness profiles that can be combined together. .. py:class:: expdisk(guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Creates an exponentially declining surface brightness profile. Inputs ------ guesses : ndarray of float Initial guesses. If a single element, then the disc profile is normalised (peak=1). If two elements then elements should be ['PeakFlux_exp','Rscale_exp']. Rscale units of arcsec. minimums : ndarray of float Minimums for the given parameters. maximums : ndarray of float Maximums for the given parameters. priors : ndarray of objects Optional- Priors for the given parameters (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameters. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = All False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: gaussian(guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Creates an gaussian surface brightness profile. Inputs ------ guesses : ndarray of float Initial guesses. If a two elements, then the gaussian profile is normalised (peak=1). If three elements then elements should be ['PeakFlux_gauss','Mean_gauss','sigma_gauss']. Units of mean/sigma are arcsec. minimums : ndarray of float Minimums for the given parameters. maximums : ndarray of float Maximums for the given parameters. priors : ndarray of objects Optional- Priors for the given parameters (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameters. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = All False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: cutoff(guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Creates an cutoff in a surface brightness profile between two radii. Inputs ------ guesses : ndarray of float Initial guesses. Two elements ['Start_cutoff','End_cutoff']. Units of arcsec. minimums : ndarray of float Minimums for the given parameters. maximums : ndarray of float Maximums for the given parameters. priors : ndarray of objects Optional- Priors for the given parameters (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameters. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = All False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:method:: eval(r, params) Evaluates a list of surface brightness profiles. Inputs ------ modellist : list of objects List of sb_prof objects, or objects that have the same methods/inputs. r: ndarray of floats Radius array, units of arcseconds. params: ndarray of floats Parameters to use in each model in the list. Returns ------- out : ndarray of floats Output combined surface brightness profile. .. py:class:: velocity_profs Circular velocity curves from physical and arbitary models that can be combined together. .. py:class:: tilted_rings(bincentroids, guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Arbitary velocity profile using the tilted ring formalism. Inputs ------ bincentroids: ndarray of float Radii at which to constrain the profile. Linearly interpolated between these points. guesses : ndarray of float Initial guesses. Size of bincentroids, units of km/s. minimums : ndarray of float Minimums for the given parameters. maximums : ndarray of float Maximums for the given parameters. priors : ndarray of objects Optional- Priors for the given parameters (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameters. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = All False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: keplarian(distance, guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Keplarian velocity profile for a point mass at r=0. Inputs ------ distance: ndarray of float Distance to the object in Mpc. guesses : ndarray of float Initial guesses. One element log10 of the central mass (in Msun). minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: mge_vcirc(surf, sigma, qobs, distance, guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Evaulate an MGE model of the potential, with or without a central point mass. Inputs ------ surf : ndarray of float Luminosity of each gaussian component, units of Lsun. sigma : ndarray of float Width of each gaussian, units of arcsec. qobs : ndarray of float Axial ratio of each gaussian. distance: ndarray of float Distance to the object in Mpc. guesses : ndarray of float Initial guesses. One or two elements [M/L and optionally log10_CentralMass in Msun]. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: arctan(guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Arctangent velocity profile. Inputs ------ guesses : ndarray of float Initial guesses. Vmax and Rturn in units of km/s and arcseconds. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: radial_barflow(guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Radial barflow from Spekkens+Sellwood. Inputs ------ guesses : ndarray of float Initial guesses. Vtangential, Vradial, BarRadius and BarPA in units of km/s, km/s, arcsec, degrees. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, args, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: sersic(distance, guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Velocity curve arising from a sersic mass distribution (spherically symmetric). Inputs ------ distance: float Distance to the object in Mpc. guesses : ndarray of float Initial guesses. Total mass, effective radius and sersic index. Units of log10(Msun), arcsec and unitless. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: mass(r, theargs, norm=1) .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, theargs, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: bulge_disc(distance, guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Velocity curve arising from a combination of a n=4 and n=1 sersic mass distributions (spherically symmetric). Inputs ------ distance: float Distance to the object in Mpc. guesses : ndarray of float Initial guesses. Total mass, effective radius and sersic index. Units of log10(Msun), arcsec and unitless. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(x, theargs, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: nfw(distance, guesses, minimums, maximums, priors=None, precisions=None, fixed=None, hubbleconst=68.0) Velocity curve arising from an NFW halo. Inputs ------ distance: float Distance to the object in Mpc. guesses : ndarray of float Initial guesses. M200 and concentration. Units of log10(Msun) and unitless. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(xs, theargs, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:class:: exponential_disc(distance, guesses, minimums, maximums, priors=None, precisions=None, fixed=None) Velocity curve arising from a razor thin exponential disc. Based on eqn 8.74 in "Dynamics and Astrophysics of Galaxies" by Bovy. Inputs ------ distance: float Distance to the object in Mpc. guesses : ndarray of float Initial guesses. Mdisk and Scaleradius. Units of log10(Msun) and arcsec. minimums : ndarray of float Minimums for the given parameter. maximums : ndarray of float Maximums for the given parameter. priors : ndarray of objects Optional- Priors for the given parameter (see GAStimator priors). precisions : ndarray of float Optional - Precision you want to reach for the given parameter. [Default = 10 percent of range] fixed: ndarray of bool Optional - Fix this parameter to the input value in guesses. [Default = False] .. py:method:: __repr__() Return repr(self). .. py:method:: __call__(xs, theargs, **kwargs) Returns the required profile. Inputs ------ x : ndarray of float Input radial array in arcseconds args : ndarray of float Input arguments to evalue the profile with .. py:method:: eval(r, params, inc=90) Evaluates a list of velocity_profs objects, returning the total circular velocity profile. Inputs ------ modellist : list of objects List of sb_prof objects, or objects that have the same methods/inputs. r: ndarray of floats Radius array, units of arcseconds. params: ndarray of floats Parameters to use in each model in the list. inc: float Inclination of system, in degrees. .. py:function:: transformClouds(inclouds, posAng=90.0, inc=0.0, cent=[0.0, 0.0], sbRad=None) Transforms skySampler cloudlets for use in KinMS. Calculate the galaxy co-ordinates of clouds from the sky plane. This MUST be used if any of the following conditions are true: inc != 0 posAng != 90 cent != [0,0] This exists as a stand-alone routine since an MCMC fit to the galaxy will likely need to run this every step, and the sampleClouds routine is computationally expensive ============ Inputs: clouds: np.ndarray The output of the sampleClouds array [x,y,I] posAng: 0=