HTML automatically generated with rman
Table of Contents
rotcur - fit parameters to velocity field using tilted rings
rotcur
in=velfie [parameters=values ...]
rotcur derives the kinematical
parameters from the observed velocity field by fitting a tilted-rings model
to the velocity field. It does a non-linear least-squares-fitting to the function:
v(x,y) = VSYS + VROT * cos(theta) * sin(INC)
- (x-XPOS) * sin(PA) -/+ (y-YPOS) * cos(PA) ?
where: cos(theta) = -----------------------------------------
r
- (x-XPOS) * cos(PA) - (y-YPOS) * sin(PA)
and: sin(theta) = -----------------------------------------
r * cos(INC)
where the radius r is measured in the plane of the galaxy:
r = sqrt( (x-XPOS)**2 + (y-YPOS)**2/cos(INC)**2 )
In the above formula v(x,y) denotes the radial velocity at rectangular
sky coordinates x and y, VSYS the systemic velocity, VROT the rotational
velocity, INC the inclination angle and theta the azimuthal distance from
the major axis in the plane of the galaxy. In itself, theta is a function
of the inclination (INC) and the position angle (PA) of the major axis.
XPOS and YPOS denote the position of the rotation center in pixels w.r.t.
0,0 being the lower left corner of the map. rotcur can fit for each ring
the 6 parameters VSYS, VROT, INC, PA, XPOS and YPOS, though any combination
of them can be fixed (see fixed=). The position angle PA of the major axis
is defined as the angle, taken in anti-clockwise direction between the north
direction on the ‘‘sky’’ and the major axis of the receding half (positive
radial velocity) of the galaxy.
Values equal the undefined value (currently
0.0) are ignored in the fit. See ccdmath(1NEMO)
or fitsccd(1NEMO)
on how
to create a velocity field with such undefined values.
By setting the keyword
fitmode= appropriately, rotcur can also be used to fit a model of pure
expansional velocities:
v(x,y) = VSYS + VROT * sin(theta) * sin(INC)
use fitmode=sin,1 in this case. currently XPOS and YPOS must be held fixed
in this case.
For a detailed discussion of this method see K. Begeman, Astr.
& Astrophys. 223, 47. (1989) and references therein.
See also velfit(1NEMO)
on the classic WWB73 (Warner, Wright and Baldwin, 1973, MNRAS, 163,163)
method.
The following parameters are recognized:
- in=image_vel
- Input velocity field map, either in image(5NEMO)
format, or table(5NEMO)
format. See also imagemode= below. No default.
- radii=r0,r1,r2,...rN
- Inner and
outer radii for N touching rings (in arcsec). Remember that the number of
radii is hence one more than the number of rings. The units= keyword can
be used to scale your physical units in the header to ‘‘human readable’’ units
(arcsec). Radii must be either all increasing or decreasing. Note that the
corresponding keyword in GIPSY has a different meaning, where they denote
the ring centers, and a different keyword WIDTHS= controls the width of
each ring. No default.
- vrot=v1,v2,...vN
- Array of rotation (or expansion) velocities,
one for each ring. If less than the number of rings is given, the trailing
ones are filled with the last supplied value.
- pa=p1,p2,...pN
- Array of position
angles (degrees), one for each ring. As is the convention, this is the
PA of the receding side of the galaxy. See also vrot=.
- inc=i1,i2,...iN
- Array
of inclinations (degrees), one for each ring. See also vrot=.
- vsys=v0
- Systemic
velocity. Only one number applies, though VSYS can vary accross rings.
- center=x0,y0
- Rotation center (grid units w.r.t. lower left corner being 0,0). Two numbers
are required. Default: center of map.
- frang=
- free angle around minor axis
(degrees), in the plane of the galaxy, from which data is excluded from
the fit (thus the total cone size of ignored data around the minor axis
is 2*frang) [Default: 20.0].
- side=
- Choose the side of the galaxy to fit the
velocity field to. Valid options are receding, approaching or both side(s).
[Default: both].
- weight=
- Choice of geometric weighting function with which
points are weighed into the least squares solution as a function of galactic
angle away from the major axis. Valid options are: uniform, cosine, and
cos-squared. [Default: cosine].
- fixed=
- List of parameters, separated by commas,
to be kept fixed during the fit. Choose any of the following: vsys, vrot,
pa, inc, xpos, ypos, although at least one parameters should be kept free.
[Default: none, i.e. all parameters free.
- ellips=
- The names of two parameters
for which to calculate an error ellips. (see fixed=). For the two parameters
it shows the major and minor axis, position angle of the one sigma deviation
ellipse. [Default: not used]
- beam=
- The beam size (FWHM, in arcsec) for beam
correction. One or two numbers required. Currently these are only used to
correct error bars for the number of independant points per beam. If not
given, each point is assumed independant. [no correction].
- dens=image_den|true|false
- Image containing the density. From this local derivatives (dN/dx)/N and
(dN/dy)/N are computed numerically, and used for an estimate of beam
smearing corrections. See also testbsc(1NEMO)
. Note that for bizarre reasons
this keyword can be used to trigger reading a weight column in tabular
input mode (Bimagemode=). [Default: required if a beam is supplied].
- wtmap=wt_image
- If given, the weights (like 1/RMS^2) of the pixels.
- tab=
- If specified, this
output table is used in append mode! This table can be piped through tabcomment(1NEMO)
and fed to ccdvel(1NEMO)
to create model velocity fields. [Default: not
used].
- resid=
- If specified, this output will either contain a residual image
(OBS-FIT) if the input file was an image, or a table with X, Y, Velocity
and Residual Velocities in each ring. Each ring will be separated from
the previous one with a simple commented line specifying which ring it
refers to. [Default: not used].
- tol=
- Tolerance for convergence of nllsqfit
[Default: 0.001].
- lab=
- Mixing parameter for nllsqfit [Default: 0.001]
- itmax=
- Maximum number of allowed nllsqfit iterations [Default: 50]
- units=
- Units
of input axes for radius and velocity. Valid options are deg, arcmin, arcsec,
rad for radius. A numeric value can also be given, in which case your image
pixel separation from the image header is multiplied by this number to
get to the ‘‘arcsec’’ that will be quoted in the tables. The units for velocity
can only be numeric, and will be the factor by which the velocities in
the map are multiplied. [Default: deg]
- blank=
- Value of the blank pixel that
needs to be ignored. [Default: 0.0].
- inherit=t|f
- Logical denoting if the initial
conditions for subsequent fitted rings should be inherited from the previous
successfully fitted ring. The fixed parameters keep of course their fixed
value. [Default: t]
- reuse=t|f
- Reuse pixels between rings. If neighboring rings
have a different geometries, it can occur that pixels will be reused. This
flag will prevent that. Note, if you don’t reuse pixels, rotcur will more
likely produce a different rotation curve if you start at outer rings,e.g.
radii=100:0:-10. [Default: t]
- fitmode=cos|sin,1
- nsigma=
- Reject outlier points
will fall outside nsigma times the dispersion away from the mean velocity
in a ring. By default, it will not reject any outliers.
- imagemode=t|f
- Image
input file mode? By default the input file is an image, alternatively a
simple ascii table with X and Y positions in columns 1 and 2, and radial
velocities in column 3, and optional errors in the radial velocity in column
4 (activated by setting dens=t). [Default: t]
- wwb73=t|f
- Use a simpler WWB73
(Warner, Wright, Baldwin 1973) linear method of fitting? [false]
The
standard output is normally not very useful; it displays, for each iteration,
the run of parameters plus the number of points and mean error in the ring.
The following awk(1)
scripts may be useful to extract information per ring:
iter. number, vsys, vrot, pa, icn, xpos, ypos, npoints, sigma_vel.
BEGIN{count=0;line="";}
{
if ($1 == "radius"){
if (count != 0){
printf("%s %s0,rad,line);
rad=$4;
}else{
count=1;
rad=$4;
}
}else{
line=$0
}
}
END{printf("%s %s0,rad,line);}
Here is an example of creating a synthetic velocity field with ccdvel,
and analysing it with rotcur:
% set r=‘nemoinp 0:100:5‘
% set v=‘nemoinp 0:100:5 | tabmath - - "100*%1/(20+%1)" all‘
% ccdvel out=map1.vel rad="$r" vrot="$v" pa=30 inc=60
% rotcur in=map1.vel radii=0:100:5 vrot=0:100:5 pa=30 inc=60 vsys=0
tab=map1.rotcur units=arcsec,1
% head map1.rotcur
radius systemic error rotation error position error inclination
error x-position error y-position error
velocity velocity angle angle
of center of center
(arcsec) (km/s) (km/s) (km/s) (km/s)(degrees) (degrees) (degrees)
(degrees) (grids w.r.t. (0,0)) (grids w.r.t. (0,0))
7.50 0.00 0.44 27.08 0.50 30.45 1.09 59.72
2.47 63.50 0.20 63.50 0.30 96
12.50 -0.00 0.24 38.54 0.28 30.66 0.41 59.27
0.94 63.50 0.13 63.50 0.20 150
...
% tail map1.rotcur
average inclination : 59.93 ( 0.043) degrees
average position angle : 30.05 ( 0.044) degrees
average systemic velocity: -0.00 ( 0.001) km/s
average x-position : 63.50 ( 0.000) grids
average y-position : 63.50 ( 0.001) grids
Failures in nllsqfit are not handled gracefully, and may error out
the program. Usage of the error= system keyword can be used to bypass such
bad rings, use with caution though and study the output.
Errorbars quoted
in the table are only an estimate since the beam size is not known. Multiply
these numbers by the square root of the number of pixels per beam to get
a more realistic estimate. If the beamsize, beam=, is given, the formulae
of Sicking (1997) was previously used to correct the errors for:
factor = sqrt(4.PI.B_x.B_y/(D_x.D_y))
but this is clearly too large, instead
factor = sqrt(PI/(4ln2).B_x.B_y/(D_x.D_y))
is now used.
Sign of the pixelsize Dx,Dy in the CCD header is ignored, and
an astronomical image is assumed. See also the reuse= keyword.
For calculations
of residuals in overlapping rings (e.g. warps) only the last ring velocity
will be used.
ccdvel(1NEMO)
, tabcomment(1NEMO)
, rotcurves(1NEMO)
,
rotcurshape(1NEMO)
, pvtrace(1NEMO)
, runvelfitss07(1NEMO)
, ccdmom(1NEMO)
,
testbsc(1NEMO)
, gal(AIPS), rocur(AIPS), velfitpatch(1NEMO)
, rotcur(5NEMO)
Begeman(1989): 1989A+A...223...47B
Warner, Wright, Baldwin (1973): 1973MNRAS.163..163W
pPXF: http://www-astro.physics.ox.ac.uk/~mxc/idl/#ppxf
2DBAT: https://github.com/seheonoh/2dbat
Bbarolo: http://editeodoro.github.io/Bbarolo/
BBAROLO: https://bbarolo.readthedocs.io/en/latest/
GALPAK3D: http://galpak3d.univ-lyon1.fr/index.html
GBKFIT: http://supercomputing.swin.edu.au/projects/gbkfit/
TiRiFiC: http://www.astron.nl/~jozsa/tirific/index.html
K. Begeman (original GIPSY Sheltran version, now also available in
C), P. J. Teuben (NEMO C version)
19/jul/83 original program KGB
9/mar/85 revision of program KGB
23/may/86 migrated to VAX-VMS KGB
27/nov/88 UNIX version pjt
8-feb-91 flushed buffers ROTCUR.DAT each write pjt
30-apr-91 moved to NEMO into C pjt
10-sep-91 documentation improved pjt
17-oct-91 added Npts to table output pjt
21-may-92 added Bob Gruendl’s rotcur awk scripts PJT
12-jun-92 added inherit=t as default PJT
13-aug-92 implemented fitmode= without XPOS,YPOS PJT
15-oct-99 compute residuals and add resid= PJT
14-mar-01 V2.5: clarifications, added nsigma= PJT
9-may-01 V2.6a: corrected error correction factor PJT
10-aug-01 clarified some differences between NEMO and GIPSY versions PJT
26-jan-02 Added unit (scale factor) for velocity too PJT
26-jun-02 V2.8: added tabular input for irregular spaced data, fixed example PJT
11-sep-02 V2.9: implemented map residual velocity field PJT
30-jan-03 V2.10: allow tables to use error in velocity PJT
2-jun-04 V2.12: finally implemented the reuse= option PJT
6-jun-20 V2.13: added wtmap= PJT
18-jan-21 V2.14: beam error factor back to standard, not Sicking PJT
Table of Contents