; FTARS
;---------------------------------------------------------------
;! Fits 1-dimensional polarization spectrum to Q/U text file
;# TASK ANALYSIS PLOT SPECTRAL ONED MODELING TV-APPL POLARIZATION
;-----------------------------------------------------------------------
;;  Copyright (C) 2020
;;  Associated Universities, Inc. Washington DC, USA.
;;
;;  This program is free software; you can redistribute it and/or
;;  modify it under the terms of the GNU General Public License as
;;  published by the Free Software Foundation; either version 2 of
;;  the License, or (at your option) any later version.
;;
;;  This program is distributed in the hope that it will be useful,
;;  but WITHOUT ANY WARRANTY; without even the implied warranty of
;;  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;;  GNU General Public License for more details.
;;
;;  You should have received a copy of the GNU General Public
;;  License along with this program; if not, write to the Free
;;  Software Foundation, Inc., 675 Massachusetts Ave, Cambridge,
;;  MA 02139, USA.
;;
;;  Correspondence concerning AIPS should be addressed as follows:
;;         Internet email: aipsmail@nrao.edu.
;;         Postal address: AIPS Project Office
;;                         National Radio Astronomy Observatory
;;                         520 Edgemont Road
;;                         Charlottesville, VA 22903-2475 USA
;-----------------------------------------------------------------------
FTARS     LLLLLLLLLLLLUUUUUUUUUUUU CCCCCCCCCCCCCCCCCCCCCCCCCCCCC
FTARS     Fits polarization spectrum to Q/U text file
INFILE
                                   Q/U spectrum input to TARS
OUTFILE
                                   RM spectrum output by TARS
NITER            0.0      4000.0   Max # fit iterations
                                      0 -> 100.
NGAUSS           1.0         4.0   Max number components
DOSPIX          -1.0         4.0   1 => fit spectral index
                                   2,3,4 => fit thickness
                                   else keep both at zero
RMSLIMIT                           RMS level above which a
                                   solution is "bad"
----------------------------------------------------------------
FTARS
Task:  Fits a single polarization spectrum for one or more components
       of total polarization, polarization angle at 0 wavelength, and
       rotation measure.  FTARS will be used on the Q and U text file
       spectrum input to TARS and the rotation measure spectrum output
       by TARS.

       The task fits up to 4 polarization components the spectrum.  It
       is based on task RMFIT, but does only a single spectrum rather
       than whole data cubes.  Just as TARS is used to test isolated
       real or model cases for FARS, FTARS is used to carry the test
       through the RM fitting stage.

       Color differentiates what is plotted.  Graphics channel 1
       (usually yellow) is used for axis labels and the data.  The
       FARS output is plotted in amplitude/phase initially.  The
       initial guess is shown as a large X in graphics channel 2
       (usually green).  This screen may be used to reset the initial
       guess.  Then the Q and U data are plotted.  Since the Q and U
       data may be on a FQID axis and so may not actually be in
       frequency order, the Q and U data are plotted as small X's.
       The initial guess is shown as a smooth curve in graphics
       channel 2 (usually green).  The fit is then shown as somewhat
       large X's in graphics channel 4 (usually cyan).

       Details of the interactive operation and options are described
       in the EXPLAIN file.  See also AIPS Memo 118, "Modeling
       spectral cubes in AIPS", May 2020 for a more detailed
       description of the use of RMFIT.

       The model fit is:
          Q(i) = SUM_j [ P_j cos ( 2 theta_j + 2 RM_j lambda(i)^2)
                  F (beta_j, lambda(i)^2) ]
          U(i) = SUM_j [ P_j sin ( 2 theta_j + 2 RM_j lambda(i)^2)
                  F (beta_j, lambda(i)^2) ]
       where P_j is the polarized flux at 1 GHz, theta_j is the
       intrinsic polarization angle, and RM_j is the rotation measure
       in radians/m/m.  F(beta_j, lambda(i)^2) is an amplitude
       function that solves for a spectral index or one of several
       possible thickness parameters.  The spectral index function
       is
          (lambda_1^2 / lambda(i)^2)^(beta_j/2)
       while the thickness functions are ('SLAB' in MODIM) and ('GAUS'
       in MODIM) and ('EXP ' in MODIM)
          sin (beta_j lambda(i)^2) / (beta_j lambda(i)^2)
          exp (-ln(2) * (beta_j lambda(i)^2 / 1.8954)^2)
          exp (-ln(2) * beta_j lambda(i)^2 / 1.8954)
       where lambda_1 is the wavelength at 1 GHz and beta_j is the
       spectral index or the thickness parameter, the latter in
       rad/m/m.  Set DOSPIX(j) to 1 for spectral index or 2,3,4 for
       SLAB, GAUS, or EXP. respectively, or set it to 0 to omit
       fitting a thickness to component j.  The 1.8954 was chosen to
       make all functions have the same half-power point.  Output
       images are written of the P_j, theta_j, RM_j, beta_j, Q0_j, and
       U0_j and their uncertainties, where
           Q0_j = P_j  cos ( 2 theta_j )
           U0_j = P_j  sin ( 2 theta_j )
       are Q and U at 1 GHz with no rotation measure.  The beta_j
       images are only written when they have been fit.

Adverbs:
  INFILE.....The input file to TARS with three columns: frequency in Hz,
             then Q in Jy/beam, then U in Jy/beam.  A semi-colon in
             column 1 means the line is a comment.
  OUTFILE....The output file from TARS containing the rotation measure
             spectrum.  The 3 columns are RM, Real part, Imaginary part
             or RM, amplitude, phase.  A semi-colon in column 1 means
             the line is a comment.  NOTE: TARS can write many RM
             spectra to a single file.  Only the first of these is used
             by FTARS.
  NITER......Maximum function evaluations during the fit of each row.
             (< 10 -> 100 for NGAUSS > 1, < 100  -> 150 for 1
             Component)
  NGAUSS.....Number of Components to fit (between 1 and 4).
  DOSPIX.....If 0.5 < DOSPIX < 1.5, solve for spectral index in the
             corresponding component, otherwise take the spectral
             index to be zero.
             If 1.5 < DOSPIX < 2.5, solve for the thickness parameter
             in the corresponding component with sin(x)/x form, or,
             if 2.5 < DOSPIX < 3.5, solve for the thickness parameter
             in the corresponding component with Gaussian form, or,
             if 3.5 < DOSPIX < 4.5, solve for the thickness parameter
             in the corresponding component with Exponential form.
             Otherwise take the thickness parameter to be zero.
             NOTE: all values of DOSPIX must be either 0 or indicate
             the same functional form.
  RMSLIMIT...If the RMS exceeds this, the solution is reported as
             "bad".  0 -> 1000.0
----------------------------------------------------------------
FTARS:  One-dimensional polarization fitting of image cubes
Documenter:  E. W. Greisen NRAO
Related Programs: IMFIT, JMFIT, SAD, SLFIT, ZEMAN, XGAUS

Initial guesses for the fit are found from the TARS output data.
When NGAUSS = 1 this can happen automatically from the row data.

For the spectrum, the first step in the fitting process is to
determine an initial guess non-interactively.  the next steps are
    1. The input FARS RM data are plotted on the TV as amplitude
       and phase.  The task initially skips step 2 and goes to step
       3 for an interactive initial guess.
    2. If either NGAUSS > 1 or a fitting has been forced, FTARS asks
       in a TV menu for instructions.  If your answer  is RE-GUESS on
       the TV, FTARS will ask you to enter a new initial guess
       (step 3).  If your answer is DO FIT, then DARS goes on to step
       4.  If the answer is QUIT, the task prints its best fit and
       exits.
    3. To enter a new initial guess, you will be prompted to position
       the TV cursor at the rotation measure of each Component (note
       that only the X position matters).  To mark each point, hit any
       button (A, B, C, or D).  FTARS then returns to step 1 above to
       plot the new guess and ask again.
    4. Once an acceptable initial guess has been found, FTARS proceeds
       to call a non-linear least squares fitting routine to determine
       the polarization parameters that appear to fit the data best.
       The answers are then checked to see if they are "reasonable" -
       negative and very large components and components with rotation
       measures outside the input data are unreasonable.
    6. The final model is added to the plot.  In the input terminal,
       you are then told if the answers are "unreasonable" and are
       shown the answers, reasonable or not.  You are then offered a
       variety of choices.  If your selected option is
          QUIT     - FTARS exits cleanly
          RE-GUESS - FTARS returns to step 1 to try again with the
                     current number of components (<= NGAUSS)
          1, 2, 3, or 4 - FTARS returns to step 1 to try again with
                     the specified number of Components (note that
                     only numbers <= NGAUSS are respected)
          HAND     - FTARS prompts you to enter the polarization
                     parameters for each component (HAND on the TV).
                     You enter in the input terminal window, the
                     total polarization in input image units, the
                     polarization position angle at zero wavelength
                     in degrees (-90 to +90), the rotation measure in
                     radians per meter squared and, if appropriate,
                     the spectral index (unitless) or RM thickness
                     (radians per meter squared).  All 3 or 4 numbers
                     must appear in one line.  It will prompt for all
                     NGAUSS components, but will change the prompt if
                     the current number of Components is < NGAUSS.
                     If you enter something besides 0 in that case,
                     then the current number of Components is
                     increased appropriately.  FTARS then returns to
                     the start of this step 6 to allow you to see if
                     you made a good guess.
          DO FIT   - After a HAND operation, a DO FIT option is
                     available.  It loops back using the hand-entered
                     parameters as the initial guess to the fitting
                     routine.
          PIXRANGE - FTARS prompts you for the Q and U display
                     ranges to re-display the data
                     Enter Qmin, Qmax, Umin, Umax in that order in
                     one line.  If you enter only two values, then
                     Umin and Umax are taken to equal Qmin and Qmax.
                     If you enter no values, then the plot reverts to
                     self-scaling.  The Q plot self-scales if Qmin >=
                     Qmax; the U plot self-scales if Umin >= Umax.
                     Note that on the re-plot, the smooth "initial
                     guess" line is now actually the result of the
                     fit.
          SHOW P&A - FTARS will switch the display between Q plus U
          SHOW Q&U   and total polarization and polarization angle.
                     The plotted "initial guess" becomes the actual
                     fit results.

