Program:      ROTCUR

Purpose:      ROTCUR derives the kinematical parameters from the observed
              velocity field by fitting tilted-rings to the velocity field.

Category:     CALCULATION, DYNAMICS, ROTATION CURVES, VELOCITY FIELDS

File:         rotcur.c

Author:       K.G. Begeman

Keywords:

   INSET=     Set (and subset) of observed velocity field.

   BOX=       Select area of velocity field [entire velocity field].

   BUNIT=     If the map units are not found in the header, they must be
              supplied by the user [KM/S].

   RADII=     Give central radii of concentric rings. Maximum number of rings
              is 512. Units are arcsec.

   WIDTHS=    Give width of rings. If number of widths is less than the
              number of radii, the last supplied width will be used for the
              rest of the rings. Units are arcsec.

   VSYS=      Give initial estimate(s) of systemic velocity(ies) in km/s.
              If the number of systemic velocities is less than the number
              of rings, the last supplied systemic velocity will be used
              for the rest of the rings.

   VROT=      Initial estimate(s) rotation velocity(ies) in km/s.
              If the number of rotation velocities is less than the number
              of rings, the last supplied rotation velocity will be used
              for the rest of the rings.

   VEXP=      Initial estimate(s) expansion velocity(ies) in km/s [0.0].
              If the number of expansion velocities is less than the number
              of rings, the last supplied expansion velocity will be used
              for the rest of the rings.

   PA=        Initial estimate(s) position angle(s) in degrees.
              If the number of position angles is less than the number
              of rings, the last supplied position angle will be used
              for the rest of the rings.

   INCL=      Initial estimate(s) inclination(s) in degrees.
              If the number of inclinations is less than the number
              of rings, the last supplied inclination will be used
              for the rest of the rings.

   CENTRE=    Initial estimates of centre of rotation in any coordinates.

   FREEANGLE= Angle around minor axis in degrees within which radial
              velocities are discarded.

   SIDE=      Which half of the velocity field should be used in the
              fitting [RECEDING and APPROACHING half].

   WEIGHT=    There are three weighting functions available: the
              UNIFORM weighting function, were all points in a
              ring have equal weights, the [COSINE] weighting function,
              were each point in a ring is weighted with |cos(theta)|,
              and the COS-SQUARED weighting function, were each point is
              weighted with cos(theta)^2.

   FIXED=     Which parameter(s) should be kept fixed [NONE]. The parameters
              are named VSYS, VROT, VEXP, PA, INCL, XPOS and YPOS. If you don't
              want to fit the inclination and the systemic velocity, you
              should type: INCL VSYS. If a fit is wanted to only one half of
              the velocity field the parameters VSYS, XPOS and YPOS are
              automatically kept fixed.

** TOLERANCE= Tolerance of least-squares fitting [0.001].

** TABLE=     Store results in a table [Y].
              The results of the fitting are saved by defaults in a GDS
              table ROTCUR**, where ** is a number from 1 to 99. If
              table ROTCUR99 is already present, table ROTCUR00 will be
              used.

   FILENAME=  Name of text file to save results [rotcur**]. Default the
              name of the file is equal to the name of the table.

Description:  This program does a least-squares-fitting to the function:

              v(x,y) = VSYS + VROT * cos(theta) * sin(INCL)
                            + VEXP * sin(theta) * sin(INCL)

                               - (x-XPOS) * sin(PA) + (y-YPOS) * cos(PA)
              with: cos(theta) = -------------------------------------------
                                                      r

                               - (x-XPOS) * cos(PA) - (y-YPOS) * sin(PA)
              and:  sin(theta) = ---------------------------------------------
                                                r * cos(INCL)

              Here  v(x,y)  denotes the radial velocity at rectangular sky
              coordinates x and y, VSYS the systemic velocity, VROT the
              rotational velocity, VEXP the expansion velocity, INCL the
              inclination angle and theta the azimuthal distance from the
              major axis in the plane of the galaxy. Theta is a function of
              the inclination (INCL) and the position angle (PA) of the
              major axis. XPOS and YPOS denote the position of the rotation
              centre.
              This program will fit for each ring the parameters VSYS, VROT,
              VEXP, INCL, PA, XPOS and YPOS.
              Note that 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 of the galaxy.

Updates:      Jul 17, 1992: KGB, Document created.
              Mar 20, 1993: KGB, fitting one half of vfield repaired.
              Mar 23, 1993: KGB, fitting of expansion velocities.
              Jan 13, 1994: KGB, keyword FILENAME= added.
              May 20, 1994: VOG, MAXRING increased from 60 to 512
              Sep 26, 1995: KGB, keyword TABLE added, modified code because
                                 of slow GDS.
              Jul 16, 1996: KGB, allow different VSYS for each ring.
              Apr 01, 1998: VOG, Increased MAXPIXEL   
              May 06, 2000: VOG, Handle rings where all pixels have same
                                 value.
              Feb 17, 2011: VOG, Removed bug with table name