================================================================================
# GUIDE TO INPUT FILE FOR GALFITM (a product of the MegaMorph project)
Including multi-band fitting, non-parametric component and MultiNest sampling.
CSL = comma separated list (must not contain any whitespace)
Where several lines are given for the same letter code, these are alternative
examples. The behaviour for multiple lines with the same letter is undefined.
================================================================================
# IMAGE and GALFIT CONTROL PARAMETER

# Input data images (CSL of FITS filenames)
# the number of input data images defines <nbands>
# the order of the bands must be maintained in all multi-band options
# the first band in the list is the 'reference band'
A) gal-r.fits,gal-g.fits,gal-i.fits

# Band labels (CSL of <nbands> labels containing no whitespace)
# (these must be unique in a case-insensitive manner)
# (can be omitted if fitting a single band)
A1) r,g,i

# Band wavelengths (CSL of values)
# (choice of wavelength units is arbitrary, as long as consistent,
#  but affects the resulting wavelength-dependence parameters)
A2) 6220,4750,7630

# Output data image block (FITS filename)
B) imgblock.fits

# Sigma image name (CSL of <nbands> FITS filenames or "none")
# (if an individual filename is specified as "none", then that sigma
#  image will be made from data; if the whole entry consists of just a
#  single "none", then all sigma images will be made from data.)
# One can also add a minimum sigma value, such that any galfit-created
# sigma image will have a minimum of that value times the sky-subtracted
# input data.
C) none
C) sig-r.fits,sig-g.fits,sig-i.fits
C) none,none,none
C) sig-r.fits,none,sig-i.fits         # perhaps unwise to do this in practice
C) none    0.1			      # min. sigma is 10% of the data flux

# Input PSF image (CSL of <nbands> FITS filenames) 
# and a single diffusion kernel (FITS filename, # or omitted)
D) psf-r.fits,psf-g.fits,psf-i.fits  #  
D) psf-r.fits,psf-g.fits,psf-i.fits  kernel.fits

# PSF fine sampling factor relative to data 
E) 1                   

# Bad pixel mask (CSL of <nbands> FITS image or ASCII coord list)
# (if an individual filename is specified as "none", then a blank
#  mask will be used; if the whole entry consists of just a single
#  "none", then all masks will be blank.)
F) none
F) mask-r.fits,mask-g.fits,mask-i.fits
F) none,none,none
F) mask-r.fits,none,mask-i.fits

# File with parameter constraints (ASCII file)
G) none              
G) constraints_filename

# Image region to fit (xmin xmax ymin ymax)
H) 1    93   1    93

# Size of the convolution box (x y)
I) 100    100

# Magnitude photometric zeropoint (CSL of <nbands> values)
J) 26.563,25.123,27.987

# Plate scale (dx dy)   [arcsec per pixel]
K) 0.038  0.038

# Display type (regular, curses, both)
O) regular             

# Options: 0=normal run; 1,2=make model/imgblock & quit
P) 0                   

# Non-parametric component
U) 0     # Standard parametric fitting
U) 1     # Turn on non-parametric component with SED homogenisation
U) -1    # Turn on non-parametric component without SED homogenisation
U) 3 0.75 25 4 40 0.0 1.0 # Customise the non-parametric schedule (n,a,b,c,d,r,t)
                           # Every n iterations, the nonparametric image is updated
                           # using a fraction of the filtered residuals, npf:
                           # npf = a * b^c / (b^c + |i - d|^c),
                           # where i is the iteration number.
                           # If n is negative do not apply SED homogenisation.
                           # r is a radius around each object centre to exclude.
                           # t is a factor to modify the filtering threshold.
   
# MultiNest
V) 0     # Use standard Levenburg-Marquardt algorithm
V) 1	 # Use MultiNest sampling algorithm (experimental and slow!)
V) 1 0 500 0.8 0.5 100000   # Customise MultiNest options:
                            # ceff,nlive,efr,tol,maxiter

# Output options
# Valid options: blank,input,model,residual -- as usual
#                component -- individual model components
#                psf -- input psf image
#                sigma -- sigma image (input or created)
#                mask -- mask image
#                nonparam -- nonparam image (if appropriate)
#                datasub -- input minus nonparam image (if appropriate)
#                nonparres -- residual minus nonparam image (if appropriate)
#                itertable -- table of parameters at each iteration
W) default  # == blank,input,model,residual and assumed if omitted
W) none     # (or any other invalid option) no images output
W) model    # or any other valid option
W) model,residual    # or any comma separated list of valid options


# THE OBJECT LIST BELOW can be however long or short as the complexity
# requires.  The user has complete freedom to mix and match the components
# by duplicating each object block.

# INITIAL FITTING PARAMETERS
#
# column 1: Parameter number
#
# column 2:
# -- Parameter 0: the allowed functions are: sersic, nuker, expdisk
#    	       	  edgedisk, devauc, king, moffat, gaussian, ferrer, psf, sky
# -- Parameter 1-10: value of the initial parameters
# -- and Parameter C0: For diskiness/boxiness (<0 = disky, >0 = boxy)
#  	     	  By default this is a CSL of the values of the parameter
#  	     	  in each of the input bands.  This may be optionally indicated
#		  by putting the word 'band' at the end of the line (before any
#  	     	  comment).
#  	     	  Only specifying a single value with multiple input bands
#  	     	  assumes the same value for all bands.
#                 This can also optionally be specified in terms of Chebyshev
#	       	  coefficients by adding the word 'cheb' at the end of the line
#  	     	  (before any comment).
#	       	  In this case one should give a CSL of at most <nbands> values
#		  corresponding to coefficients of a Chebyshev series.
# 	     	  First value of the CSL specifies the parameter value at the
#	     	  average wavelength of the input bands.
#	     	  Additional values in the CSL specify the variation in that
#	     	  parameter value with wavelength, from linear (1st-order),
#	     	  quadratic (2nd-order), up to <nbands>-order (which should be
#	     	  equivalent to fitting the value independently for each band
#		  Values omitted from the end of the CSL are assumed to be zero.
# -- Parameter Z: Outputting image options, the options are:
#              	  0 = normal, i.e. subtract final model from the data to create
#		      the residual image
#	      	  1 = Leave in the model -- do not subtract from the data
#
# column 3: This may be specified in one of two ways:
#           An integer giving the order of the Chebyshev series , e.g.,
#             0 = fixed to input value(s)
#             1 = fit a constant offset from the input value(s)
#             2 = fit a linear function of wavelength
#             3 = fit a quadratic function of wavelength, etc.
#           Note that, for >2, the input values are fit by a polynomial function
#           of the specified order before fitting begins.
#           Alternatively, one may give a CSL of at most <nbands> integers
#           indicating whether or not that coefficient is allowed to vary
#           (yes = 1, no = 0).  Values omitted from the end of the CSL are
#           assumed to be zero.
#
# column 4: comment

# Sersic function --------------------------------------------------------------

# Only this first function includes multi-band examples, but the same approach
# should work for all these functions.

 0) sersic     # Object type
 1) 300.  1    # position x [pixel]  (constant with wavelength)
 2) 350.  1    # position y [pixel]
 3) 20.0,21.0,22.0  3     # total magnitude in each band
 4) 4.30,4.40,4.5   2     # R_e in each band
 5) 5.20,5.20,5.20  1     # Sersic exponent in each band
 9) 0.30,0.30,0.30  1     # axis ratio (b/a) in each band
10) 10.0            1     # position angle (PA), same value in each band
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

 0) sersic     # Object type
 1) 300.  1    # position x [pixel]  (constant with wavelength)
 2) 350.  1    # position y [pixel]
 3) 20.00,0.2,0.0  1,1,1 cheb # total magnitude and coeffs of wavelength dependence
 4) 4.30,0.1,0.0   1,1,0 cheb # R_e and coeffs of wavelength dependence
 5) 5.20,0.0,0.0   1,0,0 cheb # Sersic exponent and coeffs of wavelength dependence
 9) 0.30,0.0,0.0   1     cheb # axis ratio (b/a) and coeffs of wavelength dependence
10) 10.0           1     cheb # position angle (PA) and coeffs of wavelength dependence
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)


# Nuker function ---------------------------------------------------------------

 0) nuker              # Object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 17.2       1       #    mu(Rb)            [surface brightness mag. at Rb]
 4) 20.5       1       #     Rb               [pixels]
 5) 1.2        1       #    alpha  (sharpness of transition)
 6) 0.5        1       #    beta   (outer powerlaw slope)
 7) 0.7        1       #    gamma  (inner powerlaw slope)
 9) 0.72       1       # axis ratio (b/a)   
10) -25.2      1       # position angle (PA)  [Degrees: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)


# Corser function ---------------------------------------------------------------

 0) corser             # Object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 17.2       1       #  mu(Rb) [surface brightness mag. at Rb]
 4) 2.5        1       #    Rb               [pixels]
 5) 0.5        1       #    alpha  (sharpness of transition)
 6) 0.5        1       #    gamma  (inner powerlaw slope)
 7) 5.0        1       #  R_e (half-light radius)   [pix]
 8) 3.5        1       #  Sersic index n (de Vaucouleurs n=4) 
 9) 0.72       1       #  axis ratio (b/a)  
10) -25.2      1       #  position angle (PA) [deg: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)


# deVaucouleur function --------------------------------------------------------

 0) devauc             # Object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 18.        1       # total magnitude    
 4) 32.        1       #     R_e              [Pixels]
 9) 0.5        1       # axis ratio (b/a)   
10) 107.       1       # position angle (PA)  [Degrees: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# Exponential function ---------------------------------------------------------

 0) expdisk            # Object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 17.0       1       # total magnitude
 4) 20.5       1       #     Rs               [Pixels]
 9) 0.3        1       # axis ratio (b/a)   
10) 25         1       # position angle (PA)  [Degrees: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# Edge-on disk function --------------------------------------------------------

 0) edgedisk            # Object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 17.0       1       # central surface brightness  [mag/arcsec^2]
 4) 10.5       1       # disk scale-height    [Pixels]
 5) 20.5       1       # disk scale-length    [Pixels]
10) 25         1       # position angle (PA)  [Degrees: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# Moffat function --------------------------------------------------------------
 
 0) moffat             # object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 16.5       1       # total magnitude     
 4) 0.5        1       #   FWHM               [Pixels]
 5) 1.5        1       # powerlaw      
 9) 0.3        1       # axis ratio (b/a)   
10) 25         1       # position angle (PA)  [Degrees: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# Ferrer function --------------------------------------------------------------

 0) ferrer             #  Component type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 10.0000     0      #  Central surface brghtness [mag/arcsec^2]
 4) 50.0000     0      #  Outer truncation radius  [pix]
 5) 4.0000      0      #  Alpha (outer truncation sharpness) 
 6) 2.0000      0      #  Beta (central slope)
 9) 0.5000      1      #  Axis ratio (b/a)  
10) 50.0000     1      #  Position angle (PA) [deg: Up=0, Left=90]
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# Gaussian function ------------------------------------------------------------

 0) gaussian           # object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 18.5       1       # total magnitude     
 4) 0.5        0       #   FWHM               [pixels]
 9) 0.3        1       # axis ratio (b/a)   
10) 25         1       # position angle (PA)  [Degrees: Up=0, Left=90]
 Z) 0                  # leave in [1] or subtract [0] this comp from data?

# The Empirical King Profile ---------------------------------------------------

 0) king               # Object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 14.9805    1       #  mu(0) 
 4) 10.1328    1       #   Rc 
 5) 51.0968    1       #   Rt 
 6) 2.0485     1       # alpha 
 9) 0.9918     1       # axis ratio (b/a)  
10) 20.7684    1       # position angle (PA) 
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# PSF fit ----------------------------------------------------------------------

 0) psf                # object type
 1) 300.       1       # position x [pixel]
 2) 357.4      1       # position y [pixel]
 3) 18.5       1       # total magnitude     
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)

# sky --------------------------------------------------------------------------

 0) sky
 1) 0.77       0       # sky background       [ADU counts]
 2) 0.000      0       # dsky/dx (sky gradient in x) 
 3) 0.000      0       # dsky/dy (sky gradient in y) 
 Z) 0                  #  Skip this model in output image?  (yes=1, no=0)


# The parameters C0, B1, B2, F1, F2, etc. listed below are hidden 
# from the user unless he/she explicitly requests them.  These can 
# be tagged on to the end of any previous components except, of 
# course, the PSF and the sky -- although galfitm won't bar you from doing 
# so, and will just ignore them.  Note that a Fourier or Bending mode 
# amplitude of exactly 0 will cause GALFIT to crash because the 
# derivative image GALFIT computes internally will be entirely 0.  If a 
# Fourier or Bending amplitude is set to 0 initially GALFIT will reset it  
# to a value of 0.01.  To prevent GALFIT from doing so, one can set it to any 
# other value.

# Wavelength-dependence can be included in all these parameters in
# the same way as described above for the 'standard' parameters.
# However, this has been less tested, so take care to check your output.

#  Bending modes
B1)  0.07      1       # Bending mode 1 (shear)
B2)  0.01      1       # Bending mode 2 (banana shape)
B3)  0.03      1       # Bending mode 3 (S-shape)

#  Azimuthal fourier modes
F1)  0.07  30.1  1  1  # Az. Fourier mode 1, amplitude and phase angle
F2)  0.01  10.5  1  1  # Az. Fourier mode 2, amplitude and phase angle
F6)  0.03  10.5  1  1  # Az. Fourier mode 6, amplitude and phase angle
F10)  0.08  20.5  1  1  # Az. Fourier mode 10, amplitude and phase angle
F20)  0.01  23.5  1  1  # Az. Fourier mode 20, amplitude and phase angle

#  Traditional Diskyness/Boxyness parameter c
C0) 0.1         0      # traditional diskyness(-)/boxyness(+)

#  PA rotation is used most often to create spiral galaxies, although
#  it can be used to fit isophotal rotations too.  Note that the parameters
#  R9 (inclination to line of sight) and R10 (sky position angle) differ 
#  very subtly (but importantly) from the classical axis ratio (q, parameter 
#  9) and position angle (PA, 10), so you need to understand the
#  distinction well.  R9 and R10 should be set to 0 if one is modeling early
#  type galaxies, because otherwise they are completely degenerate with
#  the classical q and PA.  R9 and R10 are used to project a spiral galaxy 
#  disk, assumed to be round viewed face on (i.e. inclination of 0 degrees, 
#  R9=0), to other orientations and flattening.  The disk orientation, R10 
#  is 0 when pointing up and increases counter-clockwise, just like usual.  
#  When viewed face on, the thickness of the spiral arm is controlled by 
#  classical parameter 9, i.e. the axis ratio.  When R9 and R10 are zero, 
#  everything is what you are used to with the old GALFIT.  If R9 and R10 are 
#  not 0, then the classical parameters 9 and 10 (q, PA) will be a bit 
#  un-intuitive, but perfectly internally consistent with the scheme just 
#  described.  For spiral galaxies, when and both R9 and R10 are 0, the 
#  classical position angle corresponds to the PA of the galaxy bar.  This all 
#  sounds pretty confusing because you're just reading.  But once you start 
#  playing with it it'll become more clear what I'm talking about.
#  To have GALFIT generate some example galaxies, set P=1.
#
#  Note that you can couple an arbitrary number of Fourier components with
#  coordinate rotation.  If you do so, you can create very complex,
#  multi-armed, spiral structures, and arms that have different
#  thicknesses.

R0) powerlaw	       # PA rotation function (power, log, none)
R1) 30.        1       # bar radius  [pixels]
R2) 100.       1       # 96% asymptotic radius (i.e. at 96% of tanh rotation)
R3) 275.       1       # cumul. coord. rotation out to asymp. radius [degrees] 
R4) 0.5	       1       # asymptotic spiral arm powerlaw 
R9) 0.5	       1       # inclination to L.o.S. (controls projected axis ratio)
R10) 30.       1       # sky position angle

#  The other coordinate rotation function is the log function.

R0) log                # PA rotation func. (tanh, sqrt, log, linear, none)
R1) 30.        1       # bar radius  [pixels]
R2) 100.       1       # 96% asymptotic radius (i.e. at 96% of tanh rotation)
R3) -59.9192   1       # cumul. coord. rotation out to asymp. radius [degrees]
R4) 10.	       1       # Logarithmic winding scale radius [pixels]
R9) -56.2757   1       # Inclination to L.o.S. [degrees]
R10) 157.9937  1       # Sky position angle

#  Create a truncation by multiplying a profile function with a
#  hyperbolic tangent transition.  The break radius, r_break, is defined
#  as the radius which has 99% the flux of an original function, whereas
#  the softening radius, r_soft, is where the function has only 1% of the 
#  flux of the original.  Under normal circumstances, r_break < r_soft.  
#  On the other hand, r_break > r_soft is also possible.  This happens,
#  for instance, if one wants to generate a ring model with a hyperbolic
#  tangent truncation in the inner region.
#
#  When multiple profiles are linked, mathematically they are linked by:
#    f_net(r) = Sum{f_inner(r)} * (1-s(r)) + Sum{f_outer(r)} * s(r).
#  The transition radii are shared by the inner & outer functions in the
#  sense that r_break is where inner profiles reaches 99% of their normal 
#  flux; it is also the radius where the outer functions reaches 1% of their 
#  own normal flux, and vice versa.  When this happens, be careful that the 
#  sense of the break vs. asymptotic radius are in the correct sense, or 
#  else things can get very confusing very quickly.  So as to minimize 
#  confusion when linking, consistently use r_break < r_soft.
#  
#  When linking two profiles, the flux parameter for the inner component 
#  is the surface brightness at the original (untruncated) effective 
#  radius "Re" for the Sersic profile instead of the total flux, or central
#  surface brightness for all other profiles (gaussian, moffat, etc.).
#  The outer component flux normalization is the surface brightness
#  at the **break** radius.  
#
#  Note that the Fourier and bending modes can operate on the truncation
#  parameters independently of the light component that the truncation 
#  parameters are modifying.
#  
#        
#        v---  r_break ---v
#        
#        ____           ___
#       /                  \
#   ___/                    \___
#
#      ^----   r_soft   ----^
#
#  There are 4 kinds of truncation modes, designated by, "radial", "radial-b",
#  for Type 1, versus "radial2", and "radial2-b" for Type 2.  
#
#  The difference between Type 1 and Type 2 is in the definition of 
#  parameters T3 and T4.  For Type 1, T3 is "Break radius", and T4 is
#  the "Softening length", i.e. (R_break) and (Delta R).  For Type 2, T3 is 
#  still the "Break radius", but T4 is the "Softening radius", i.e. 
#  inner and outer break radius.  
#
#  As for Type "b" versus Type "a" (i.e. non-b), the difference is that
#  "Type a" is intended for spiral structures, i.e. the truncation shape 
#  (axis ratio and PA) are tilted and rotated by the same angles as the 
#  spiral arms.  Whereas for "Type b" truncations, the shape parameters
#  refer to how they appear in projection, i.e. in the plane of the sky.
#  For non-spiral models, there's no difference between "Type a" and
#  "Type b".

# In principal wavelength-dependence can be included in all the
# parameters of these truncation components in the same way as
# described above for the 'standard' parameters.  However, this is
# currently not tested and should be considered unsupported.

# Truncated Sersic function
 0) sersic             # Object type
 1) 300.  350.  1 1    # position x, y        [pixel]
 3) 20.00      1       # total magnitude    
 4) 4.30       1       #     R_e              [Pixels]
 5) 5.20       1       # Sersic exponent (deVauc=4, expdisk=1)  
 9) 0.30       1       # axis ratio (b/a)   
10) 10.0       1       # position angle (PA)  [Degrees: Up=0, Left=90]
Ti) 5                  # Inner truncation by component 5
To) 2                  # Outer truncation by component 2
 Z) 0                  # leave in [1] or subtract [0] this comp from data?

# Object number: 5
T0) radial                 #  truncation 
T1) 200.  150.  1 1        #  Centroid of truncation function (optional)
T4) 4.4179      1          #  Break radius (99% normal flux)    [pixels]
T5) 9.1777      1          #  Softening length (1% normal flux) [pixels]
T9)  0.7        1          #  Axis ratio (optional)
T10)  -32.      1          #  Position angle (optional)
F1) 0.1  30     1          #  Fourier mode1 (now modifying truncation)