spider scan calibration routines

Last modified: Thu Mar 16 13:07:12 2023.


List of Routines


Routine Descriptions

ALLBEAMS_EVAL - RETURN MAPS OF ALL BEAMS

[Next Routine] [List of Routines]
NAME:
allbeams_eval -  Return maps of all beams
PURPOSE: Return maps of all beams (the Stokes I beam, the sidelobe, the
         squintbeam, the squashgeam0 on a 60 by 60 grid.

CALLING SEQUENCE:
	ALLBEAMS_EVAL, nterms, arcmin, pixelsize, azarray, zaarray, $
   			stripfit, b2dfit, $
   			mainbeam, sidelobe, $
   			squintbeam, squashbeam, totalbeam, $
  			    mainbeam_integral, sidelobe_integral, $
               squintbeam_integral, squashbeam_integral, totalbeam_integral

INPUTS:
   NTERMS: the nr of terms in the Fourier expansion of the first
			sidelobe to include. We set nterms=6, which allows it to get the
			triangle contribution.

   ARCMIN: the nominal HPBW in arcmin used in the observing
			pattern. From HDR1INFO[28, *].

PIXELSIZE: the size of the pixels. From MAKE_AZZA

  AZARRAY: ZAARRAY, the az,za array used for observing and for
			which the beam is calculated. From MAKE_AZZA

 STRIPFIT: the STRIPFIT array for this particular pattern. See
			BEAM_DESCRIBE 

   B2DFIT: the main beam description array for this particular
			pattern. See BEAM2D_DESCRIBE.

OUTPUTS:
   MAINBEAM:The Stokes I main beam (no sidelobes). This array is
			 fltarr(60,60,4), giving the possibility of 4 stokes parametres,
            but only Stokes I is represented.

   SIDELOBE[60,60,4]:The sidelobe in all four stokes parameters

   SQUINTBEAM[60,60,4]: The squint beam in the latter 3 Stokes
						parameters, so there are only three squintbeams.

   SQUASHBEAM[60,60,4]:The squash beam in the latter 3 Stokes
				parameters, so there are only three squashbeams.

   TOTALBEAM: SQUINTBEAM+ SQUASHBEAM+ SIDELOBE

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/allbeams_eval.pro)


ALLRESPONSE_EVAL - EVAL RESPONSE ALL BMS TO BRIGHTTEMP GRAD.

[Previous Routine] [Next Routine] [List of Routines]
NAME:
allresponse_eval - eval response all bms to brighttemp grad.

 PURPOSE: Evaluate the response of all beams to a brightness temperature
temperature gradient and second derivative. Evaluate the angle to which
the beams are most sensitive. 

 CALLING SEQUENCE:
ALLRESPONSE_EVAL, arcmin, stripfit, b2dfit, $
   mainampl, mainangl, sideampl, sideangl, $
   squintampl, squintangl, squashampl, squashangl, $
   totalampl, totalangl
   
INPUTS:
   ARCMIN, the nominal HPBW in arcmin used during observing. From
HDR1INFO[28, *].

       STRIPFIT, the STRIPFIT array for this particular pattern. See
BEAM_DESCRIBE

       B2DFIT, the main beam description array for this particular
pattern. See BEAM2D_DESCRIBE.

OUTPUTS:
   MAINAMPL[ 2] Amplitude of main beam response to the
two derivatives.

   MAINANGL[ 2] The position angle of two derivatives at which the
main beam has max response

   SIDEAMPL[ 4,2] and SIDEANGL[ 4,2] Amplitude and position angle
for the sidelobes, all four Stokes parameters,

   SQUINTAMPL, SQUINTANGL, SQUASHAMPL, SQUASHANGL[ 4,2], amplitude
and position angle of the squint and squash beams, the latter 3 Stokes
parameters only.

   TOTALAMPL, TOTALANGL: the totalbeam (=sidelobe+ squint+ squash)
responses to the polarization.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/allresponse_eval.pro)


AVG_SRCPHASE - CALC THE VECTOR AVERAGE ANGLE OF A SET OF INPUT ANGLES

[Previous Routine] [Next Routine] [List of Routines]
NAME:
 avg_srcphase - calc the vector average angle of a set of input angles

CALLING SEQUENCE:

  avg_srcphase, sourcephase, srcphase

INPUTS:

   SOURCEPHASE: the array of phases to average, units ***RADIANS***

OUTPUTS:

      SRCPHASE: the average phase, units ***RADIANS***.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/avg_srcphase.pro)


BEAM1DFIT - 1D FITS TO STOKES I,Q,U,V FOR EACH STRIP

[Previous Routine] [Next Routine] [List of Routines]
NAME:
beam1dfit - 1d fits to stokes i,q,u,v for each strip
SYNTAX: beam1dfit,nrc, beamin_arr, beamout_arr,stkoffsets_chnl_arr, $
                  byChnl=byChnl

DESCRIPTION
Fit the total intensity (the one-dimensional fit) separately for each
strip using ifit_allcal.pro, which fits three gaussians--one for each
sidelobe and one for the main beam. The main beam gaussian has a
skewness parameter.

Then, using the results from the x+y fit,  fit the 1-d beam for all
polarized stokes params for all strips using sqfit_allcal.pro, which
fixes all parameters for the sidelobe gaussians except for the height
and for the main beam gaussian includes squint and squash terms.

Note that no input angles other than the total offset from center 
are required to do these fits, so the only angular variable is
totoffst[ ptsperstip,4]

INPUTS: 
       NRC: the pattern number.
    BEAMIN: the input data structure, from which we extract...
HPBW_GUESS: the value of the HPBW to use as the initial guess
            in the nonlinear 1d Gaussian fits to each strip cut across 
            the beam.
TOTOFFSET[ ptsperstrip, nrstrips]: the total angular offset from
            center along each strip. UNITS ARE ARCMIN, in contrast to the
            original version of this program in which units were the assumed
            hpbw. nrstrips is 4 in original usage.

STOKESOFFSET[4,ptsperstrip,nrstrips]: the continuum Stokes parameters [4] 
            for each point in each of the nrstrips strips.

KEYWORDS:

    CHNLS: if set, does each chnl, reading inputs from
           beamin_arr.stkoffsets and writing to beamout_arr.stripfit_chnls. 
           Otherise it does the continuum, which are the average over many 
           chls. 
OUTPUTS:

BEAMOUT_ARR: INTO WHICH WE INSERT EITHER STRIPFIT (FOR THE CONTINUUM)
             OR STRIPFIT_CHNLS( VALUES FOR EACH AND EVERY SINGLE CHANNEL)...

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/beam1dfit.pro)


BEAM2DFIT - 2D LEAST SQUARES FIT FOR THE BEAM PLUS SIDELOBES.

[Previous Routine] [Next Routine] [List of Routines]
NAME:
beam2dfit - 2d least squares fit for the beam plus sidelobes.

INPUTS:

   HPBW_GUESS, the hpbw guessed value for ls fit.  all internal offset
angles, and some outputs, are specified in units of hpbw_guess. 

   
   BEAMIN, the structure containing input data.  Inputs from
structure BEAMIN that are used include... 

   BEAMIN.AZENCODERS, the set of [60,4] az encoder readings for the pattern.
   BEAMIN.ZAENCODERS, the set of [60,4] za encoder readings for the pattern.


   BEAMIN.AZOFFSETS, BEAMIN.ZAOFFSETS, the set of [60, 4] az, za
offsets in the pattern.  UNITS ARE ARCMIN, NOT HPBW; THIS IS A CHANGE
FROM ORIGINAL VERSION.  (the array size [60, 4] is illustrative; sizes
are determined from the input data.) these are immediately copied to the
(azoffset, zaoffset) arrays, which are in units of hpbw_guess, and which
are used in the calculations. 


   TOTOFFSETS, the set of [60, 40] total offsets (same for each
strip...runs negative to positive).  UNITS ARE ARCMIN, NOT HPBW; THIS IS
A CHANGE FROM ORIGINAL VERSION.  (the array size [60, 4] is
illustrative; sizes are determined from the input data.)

   STOKESOFFSET_CONT, the total system temp (antenna plus rcvr) for
each pattern point.  stokesoffset_cont[ 4, 60, 4]; again, sizes are
determined from the input data. 
   
   STRIPFIT, the array of outputs from least square fitting each
strip previously in BEAM1DFIT.

KEYWORD: 
   SQUOOSH. set it to include squoosh in the fit.

OUTPUTS

   B2DFIT is the output array containing the solved-for ls
coefficients, size [50,2]. Not all the elements are used to reserve room
for future expansion if required (God forbid!). In the second element,
The first element is the value, the second its formal error from the
fit.

   b2dfit = fltarr[ 50,2]      ;NORMALLY 0 IS VALUE, 1 IS ITS ERROR
   b2dfit[ 0,*] = [ tsys,      ;OFFSRC TSYS, KELVINS
   b2dfit[ 1,*] = [dtsys_dza,  ;DTSYS/DZA, KELVINS/DEG
   b2dfit[ 2,*] = [ tsrc,      ;SRC DEFLN, KELVINS
   b2dfit[ 3,*] = [ az_bmcntr, ;AZ OFFSET OF SCAN CENTER
                       FROM TRUE BEAM CENTER, ARCMIN
   b2dfit[ 4,*] = [ za_bmcntr,     ;ZA OFFSET, AS FOR AZ
   b2dfit[ 5,*] = [ bmwid_0,   ;AVG HPBW, ARCMIN  ***SEE NOTE BELOW***
   b2dfit[ 6,*] = [ bmwid_1,   ;(MAXHPBW-MINHPBW)/2, ARCMIN *SEE NOTE*
   b2dfit[ 7,*] = [ phi_bm, ]  ;PA OF HPBW MAJOR AXIS, DEGREES
   b2dfit[ 8,*] = [ alpha_coma,    ;COMA AMPLITUDE, UNITS OF HPBW
   b2dfit[ 9,*] = [ phi_coma,  ;PA OF COMA LOBE, DEGREES 
   b2dfit[ 10,0] = [ sigma,    ;SIGMA OF DATAPOINTS, KELVINS
   b2dfit[ 10,1] = hpbw_guess]     ;nominal HPBW ARCMIN USED IN
       OBSERVING, and the guess for the nonlinear Gauss fit

   b2dfit[ 11,*] = [ beam_integral,;MAIN BEAM INTEGRAL, K ARCMIN^2

THE FOLLOWING ARE DEFINED IN PLOT2D_BEAMFITS.PRO
   b2dfit[ 12,*] = adopted source flux. THIS IS SET EQUAL TO 1 JY
       UNLESS IT IS SPECIFIED IN ADVANCE.
   b2dfit[ 13,*] = hgt[0]: the zeroth fourier coeff of teh first sidelobe
       (same as the average sidelobe level) UNITS ARE TSRC, i.e. 
       this is ratio of sidelobe hgt to mainbeam peak 
   b2dfit[ 14,*] = eta_mainbeam for the adopted source flux
   b2dfit[ 15,*] = eta_sidelobe for the adopted source flux
   b2dfit[ 16,*] = kperjy for the adopted source flux
   b2dfit[ 17,*] = the center freq in MHz (set OUTSIDE this program)

THE PREVIOUS ARE DEFINED IN BEAM2D_PLOT.PRO

   .   .   .   .

   b2dfit[ 18,0] = the mean PA for the pattern
   b2dfit[ 18,1] = ptsperstrip, nr points per strip
   b2dfit[ 19,*] = [meanaz, meanza], meam az and za for the pattern
obtained by averaging all az's and za's for the 240 observations in the
pattern. 


   For nstk 1 to 3 (Stokes QUV), indices 20 to 29, 30 to 39, 40 to
49 we have the 3 polarized Stokes parameters, each of which is defined
as follows:
   b2dfit[ 10+ 10*nstk+ 0,*] ;ZERO OFFSET, KELVINS
   b2dfit[ 10+ 10*nstk+ 1,*]  ;Doffset/DZA, KELVINS/DEG
   b2dfit[ 10+ 10*nstk+ 2,*]  ;SRC DEFLN, KELVINS
   b2dfit[ 10+ 10*nstk+ 3,*]  ;SQUINT MAGNITUDE, ARCMIN
   b2dfit[ 10+ 10*nstk+ 4,*]  ;SQUINT PA (IN AZ/ZA SYSTEM), DEGREES
   b2dfit[ 10+ 10*nstk+ 5,*]  ;SQUASH MAGNITUDE, ARCMIN 
   b2dfit[ 10+ 10*nstk+ 6,*]  ;SQUASH PA, DEGREES
   b2dfit[ 10+ 10*nstk+ 7,*]  ;SQUASHAVG (0th TERM IN 'FOURIER EXPANSION'
                   OF SQUASH, ARCMIN 

***********IMPORTANT NOTE ABOUT DEFINITION OF POSITION ANGLE PA******

   POSITION ANGLES (PA'S) ARE IN THE AZ, ZA COORDINATE SYSTEM AND
ARE MEASURED FROM AZ=0 TOWARDS POSITIVE ZA.
THIS IS DEFINITELY NOT THE CONVENTIONAL ASTRONOMY DEFINITION OF PA!!
RATHER, IT CORRESPONDS TO PHI IN THE AOTM WRITEUP.

************** IMPORTANT NOTE ABOUT DEFINITION OF WIDTH **************
   THE INPUTS AND OUTPUTS TO THIS PROC FOR THE WIDTH ARE HPBW,
NOT 1/E. THUS B2DFIT IS HPBW IN ARCMIN.
   HOWEVER, WITHIN THIS ROUTING, THE WIDTHS ARE CONVERTED TO
1/E.

MODIFICATION HISTORY:
   Definition of alpha_coma changed 17 oct 2000. see
   Bug in hpbw conversion fixed 29 oct 2000. affects defn of alpha_coma.
   calculate meanaz using vector avg, 30 oct 2000.
   squash_avg term added (vectors from s2wdfit have 8 els instead
   of 7; include results in b2dfit). 31jul03

g2dcurv_allcal.pro, mainbeam_eval.pro

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/beam2dfit.pro)


BEAMOUT_AZZAPA

[Previous Routine] [Next Routine] [List of Routines]
NAME:
 beamout_azzapa
PURPOSE:

fill beamout structure with az, za, and beamout values for center of
scans... 

INPUTS: from beamin:
   beamin.azencoders, 240 az encoder values
   beamin.zaencoders, 240 za encoder values

OUTPUTS:
   az, za, and pa for each scan center in the structure beamout. 
all values in degrees. 

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/beamout_azzapa.pro)


CALC_BEAM2D- EVALUATE PROPERTIES OF THE 2D BEAM AND LOAD B2DFIT.

[Previous Routine] [Next Routine] [List of Routines]
NAME:
calc_beam2d- evaluate properties of the 2d beam and load b2dfit.
INPUTS:

   from structure BEAMIN

   STRIPFIT, the 1-d stripfit ls fit coefficients from BEAM_DESCRIBE

   BEAMOUT, a structure which contains...
       B2DFIT, the 2-d ls fit coefficients from BEAM2D_DESCRIBE.
       fhgt, fcen, fhpbw, generated in this proc.

   SOURCEFLUX, used for getting kperjy

OUTPUTS: ALL CALCULATED QUANTITIES OF SIGNIFICANCE ARE PUT INTO B2DFIT.

   SIDELOBE_INTEGRAL, the integral of the sidelobe pattern in units
of Kelvins arcmin^2.

   MAINBEAM_INTEGRAL, the integral of the mainbeam pattern in units
of Kelvins arcmin^2.

   SIDELOBE, the sidelobe pattern on the original 60 X 60 observing
grid.

   MAINBEAM, the mainbeam pattern on the original 60 X 60 observing
grid. 

   KPERJY, the kperjy of the source. here Kelvins is Stokes I/2.

   FHGT, the ratio of the heights of first sidelobe to main beam.

   ETA_MAINBEAM, the main beam efficiency

   ETA_SIDELOBE, the sidelobe efficiency
   
KEYWORDS:

   NTERMS, the number of terms to use in reconstructing the
sidelobe structure from the Fourier fit to the sidelobe structure. Use
of NTERMS=8 forces the reconstruction to duplicate the measured points
because there were 8 measured data points (2 on each end of the strip,
so they are located every 45 deg around the circle). Use of NTERMS<8 is
equivalent to least squares fitting those 8 points with a Fourier series
having the chosen value of NTERMS, which means that the fitted curve
will not go through the points. Seeing as the 'measured points' are really
the result of model fits to the sidelobes with Gaussians, and that these
fitted parameters are not perfect, I don't believe using NTERMS=8 is
justified. Rather, some smaller value is better. Currently I am using 6.
This needs some investigation by looking at the sidelobe coefficients,
the original strip data, and the reconstructed results using different
values of NTERMS to see what's appropriate. 

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/calc_beam2d.pro)


CROSS3_NEWCAL - INPUT, CALIBRATE, AND 1D PROCES A PATTERN

[Previous Routine] [Next Routine] [List of Routines]
NAME:
cross3_newcal - input, calibrate, and 1d proces a pattern
SL is a scanlist of the entire file. slPatInd is and index into
sl[slpatInd] that is the scan for the start of the pattern we want to
process.

Processs a single pattern. Input, phase and ;intensity calibrate. 
PROCESSES ALL 242 SCANS IN THE PATTERN, CALIBRATING THE LATTER 240 IN
TERMS OF THE FIRST TWO (WHICH ARE CALON AND CALOFF). THE CAL IS USED FOR
TWO PURPOSES:

   1. AMPLITUDE CALIBRATION. This assumes that the gain of each
spectral channel is proportional to the power in the cal-off spectrum.
It then determines the temperature scale from the cal deflection and the
values of tcalxx, tcalyy, which are the cal values as found either from
phil's file or inputted by the user.

   If the cal value in the arrays tcalxx_board and tcalyy_board is
negative, it uses the value from Phil's files. It it is positive, it
uses that value. 

   Special case: if the rcvr is 430ch, then it uses tcalxx_430ch
and tcalyy_430ch.

   2. PHASE CALIBRATION. If the system has a correlated cal
(nocal=0 as determined from the routine getrcvr.pro), then it assumes
the cal phase is zero and corrects all 240 phases to the cal phase. this
phase calibration includes a linear variation with frequency, using
carl's patented routine! 

   If the system does not have a correlated cal, then it fits the
slope of the phase with frequency and corrects all 240 observations for
that, but does not fit the zero of phase--it calculates the zero point
phases of the 240 points from their xy and xy data without subtracting
any other phase (as opposed to the correlated cal case, for which it
subtracts the zero point of the cal phase).

INPUTS:

   There are several inputs passed through the structure SCNDATA.  
These are of the nature of which points in the scan to use for on source 
and for off source values, and which points have the cal on and off. 
Here are some of the more important:

Variables contained in the structure SCNDATA:
   INDXCALON, INDXCALOFF: the index numbers within the pattern for
calon and caloff, equal to [0,0] and [1,1]. These are arrays because in
principle there could be more than one index used for calibration; if
there is only one, avoiding errors requires making them arrays with two
identical elements.

   ONSCANS, OFFSCANS: Arrays of index numbers within the pattern
used for calculating source deflections, onsource - offsource.

   DPDF, the initial guess for the phase slope in radians per MHz.
The slope arises primarily in the i.f. and its sign depends on whether
the final sideband is upper or lower; the program determines this from
whether the spectrum is flipped. In the unlikely event that the slope
becomes determined from other causes, this may need to be
changed--however, the fitting routine is pretty robust and even a guess
with the wrong sign isn't very serious in most (all?) cases.

   TCALXX_BOARD, TCALYY_BOARD: See the above discussion of
AMPLITUDE CALIBRATION.

   TCALXX_430CH, TCALYY_430CH: See the above discussion of
AMPLITUDE CALIBRATION.

   CHNLS, PHASECHNLS, GAINCHNLS: the array of channels over which
sums are taken to determine continuum power. Normally equals
indgen(128). I always set these arrays equal to each other; some
calculations use different ones, which allows different channels to be
used for phase and intensity sums, but i'm not sure if this is
consistently defined. If you want to make these arrays different, you
need to check to see how the program differentiates between them, and
whether it does so consistently. 

 POLBAD: int   ==0 --> polA bad, ==1 --> polB bad. If bad then copy other
                  other pol to this pol. This makes Q 0 and u,v junk
 missingcal: int  if set then this receiver has no cal. normalize to
                  Tsys.. you can still get pointing and sefd
excludeBrds[n]  :  do not process these boards. brds are 0..3
                   -1 in first entry ignores keyword

   
OUTPUTS:

   beamin_arr[curFitInd:curFitInd+numBrdsProcessed-1]
         This holds the calibrated data and the calibration info
   hb_arr[4,curFitInd+numBrdsProcessed-1]:{hdr} 
         Holds the header for the first sample of each strip.
   numBrdsProcessed: long  number of boards processed on this call.
                       You should increment curFitInd by this amount.

   RETURNSTATUS: status of reading the corfile. 1 if normal; 0 if
               eof; -1 if some other problem

KEYWORDS:

   set CUMCORR to excise interference from the calon/caloff spectra

   TOTALQUIET suppresses all printed output for duncan

HISTORY:
   scannr0_first = scannr0 inserted in monitor mode 6 jun 01. 
   comment inserted just above :
   09nov02 - pjp.. changed to use scanlist processing.
   01jan05 - pjp.. passed in tcalxx_430ch, yy_430ch
   13aug06 - pjp.. added han keyword in case strong rfi.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/cross3_newcal.pro)


DFT - DIRECT FT.

[Previous Routine] [Next Routine] [List of Routines]
NAME:
dft - DIRECT FT. 
	USES IDL CONVENTION FOR FORWARD/REVERSE WHICH IS OPPOSITE NM'S

CALLING SEQUENCE: dft, inputx, inputy, outputx, outputy, [/inverse]

INPUTS:
	inputx, x-value at which input sig is sampled
	inputy, amplitude of input sig
	outputx, x-values at which you want the ft computed.
KEYWORD: 
	INVERSE. if set, does the idl inverse equivalent.
OUTPUTS
	outputy, ft of sig evaluated at outputx. complex.

TO TEST AND COMPARE WTIH THE FFT: see the attached 

	pro dfttst_2


IN PARTICULAR, NOTE HOW THE ORDER OF THE INPUT AND OUTPUT ARRAYS FOR THE
FFT MUST BE REARRANGED SO AS TO BE MONOTINICALLY INCREASING

MODS
	11nov04, carlh added check to see which method is faster in for loop.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/dft.pro)


FFTEVAL-

[Previous Routine] [Next Routine] [List of Routines]
NAME:
ffteval-
PURPOSE: evaluate an inverse Fourier series using a smaller number freqs
than the original FT calculated. For example, we began with 8 points
around a circle, called pts, and got FHGT= fft( pts), there are 8 FHGT
values at 8 different frequencies. the original data were real so the 8
FHGT values are hermitian. now we can calculate the original points
using all the 8 frequencies or just some of them, e.g. 6, 4, or 2 (the
number must be even). NTERMS is the number of frequencies, equal to 8 to
use all the original freqs, or equal tl 6, 4, or 2 to retain only the
low freq terms.

CALLING SEQUENCE:
	RESULT= FFTEVAL( nterms, fhgt, theta)

INPUTS: 
	NTERMS is nr of frequencies actually employ. NTERMS must be even
and less than or equal to the number of terems in fhgt.

	FHGT, the complex array of Fourier coefficients. must have 2^N elements.

	THETA, the angle RADIANS at which to evaluate the Fourier series.

OUTPUT: 
	VALUE of the Fourier series at angle THETA.  As originally
calculataed it's complex, but the imag part is zero so it is conferted
to float. 

HISTORY: Written by Carl Heiles, original version totally wrong unless
nterms was equal to n_elements( fhgt). redone 31 may 2005. redone
again 1 jun. redone again 2 jun. redone finally (we hope!) 6 mun.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/ffteval.pro)


FT_SIDELOBES_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
ft_sidelobes_newcal
PURPOSE: generate the fourier coefficients that describe the hgt, cen,
and wid of the Gaussians that describe the sidelobes.

CALLING SEQUENCE:
FT_SIDELOBES, stripfit, b2dfit, fhgt, fcen, fhpbw, nstk=nstk

INPUTS:

   STRIPFIT, the results of the 1-d strip fits described and
derived in BEAM_DESCRIBE.

   B2DFIT, the results of the 2-d pattern fits described and
derived in BEAM2D_DESCRIBE. 

OUTPUTS:

   FHGT, FCEN, FHPBW: the Fourier coefficients of the sidelobe hgt,
cen, and width. These are complex and have dimension equal to twice
the number of strips--e.g. 4 strips means they have dimension 8 (one 
for each end). .

KEYWORD:
   NSTK: IF SPECIFIED, RETURNS THE BEAM FOR THE PARTICULAR STOKES
   PARAMETER. 0,1,2,OR 3. IF NOT specified, it returns Stokes I.

******************* IMPORTANT NOTE NR 1 **************************
   the angle units are ***ARCMIN***, so when the FT coefficients
are used to obtain the sidelobe the angular units must be ARCMIN.

******************* IMPORTANT NOTE NR 2 **************************
   the hgt, cen, and hpbw quantities must be ordered properly so
that as you go around the circle they increase sequentially. The
ordering below is correct for the 4-strip Arecibo pattern. 

   This program will also do a 2-strip pattern. However, you must
make sure that the ordering of the points is correct!
******************************************************************

COMMENT: all widths are HPBW

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/ft_sidelobes_newcal.pro)


G2DCURV

[Previous Routine] [Next Routine] [List of Routines]
NAME:
   G2DCURV

PURPOSE:
    Calculate multiple (N) Gaussians + offset + slope (wrt ZAencoder).

CALLING SEQUENCE:
    

INPUTS:
FIRST, THE OBSERVED DATA FROM THE TELESCOPE:
   DELT_ZAENCODER is the zenith angle from the encoder MINUS an offset
       which should be the mean of the za_encoder values 
       (makes the fit meaningful)
   AZ_SCAN is the set of az offsets from assumed beam center
   ZA_SCAN is the set of za offset from assumed beam center 

NEXT, THE PARAMETERS OF THIS FUNCTION:
FOR THE CENTRAL GAUSSIAN:
   TSRC is the height of the beam center
   AZ_BMCNTR is the az offset of assumed beam cntr from true beam cntr
   ZA_BMCNTR is the za offset of assumed beam cntr from true beam cntr
   BMWID_0 is the constant term in the beamwidth
   BMWID_1 is the coefficient of the 2-phi term in the beamwidth
   PHI_BM is the angle of the major axis wrt az=0 of the ellipse 
       of the beam. ***RADIANS***
   ALPHA_COMA is the magnitude of the coma lobe
   PHI_COMA is the angle of the coma lobe wrt az=0. RADIANS.
   TSYS is the off-source system temp
   DTSYS/dza is the derivative of TSYS wrt ZA

FOR THE SIDELOBES:
   STRIPFIT contains the amplitudes and other parameters of the
sidelobes--see comment at 'tout =...' below. 

OUTPUTS:
   TOUT, the array of output temps.

NOTE REGARDING COMA: We assume the coma lobe to be parameterized by a
DIRECTION, called PHI_COMA, and a MAGNITUDE, called ALPHA. 

RELATED PROCEDURES:
   
HISTORY:
   Written by Carl Heiles. 23sep00.
   17 oct: modified sign of alpha and what it is divided by.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/g2dcurv_allcal.pro)


G2DFIT

[Previous Routine] [Next Routine] [List of Routines]
NAME:
   G2DFIT

PURPOSE:

   A 2d nonlinear least squares fit to central beam including coma
and elliptical beam with arbitrary position angle. Called from
BEAM2D_DESCRIBE.

   ***IN CONTRAST TO MOST OF OUR GFITTING ROUTINES, HERE THE
INPUTS AND OUTPUTS ARE 1/E WIDTHS!!!***

INPUTS:

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/g2dfit_allcal.pro)


GCURV_ALLCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
gcurv_allcaL
PURPOSE:
    Calculate the function used in beam and sidelobe fitting, which
    consists of 3 gaussians plus zero plus slope.
    The FIRST Gaussian is for the main beam and includes a coma parameter, alpha.

CALLING SEQUENCE:
    GCURV_ALLCAL, xdata, tsys1, slope1, alpha1, hgt1, cen1, wid1, tfit

INPUTS:
     xdata: the x-values at which the data points exist.
     tsys1: the estimated constant zero offset of the data points.
     slope1: the slope.
     hgt1: the array of N estimated heights of the Gaussians.
     cen1: the array of N estimated centers of the Gaussians.
     wid1: the array of N estimated halfwidths of the Gaussians.

OUTPUTS:
     tfit: the calculated points.


HISTORY:
   Written by Carl Heiles. sep00

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/gcurv_allcal.pro)


GET_OFFSETS_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
get_offsets_newcal
GET THE ANGULAR OFFSETS, SPECTRA, AND SPECTRUM-INTEGRATED POWERS 
   FOR THE FOUR STRIPS IN ARRAYS OF [*, 4]
THE UNITS OF ANGULAR OFFSETS ARE ---> NO LONGER <--- THE SPECIFIED HPBW.
THE UNITS OF ANGULAR OFFSETS ARE ---> ARCMIN <--- .

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/get_offsets_newcal.pro)


GFIT_ALLCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
   GFIT_ALLCAL

PURPOSE:

    FITS 1, 2, OR 3 GAUSSIANS. THE FIRST GAUSSIAN IS MEAND TO REPRESENT
THE MAIN BEAM AND HAS A SKEWNESS PARAMETER, ALPHA. THE OTHER ARE
STANDARD GAUSSIANS WITH HGT, CEN, WID. ALSO, THE ZERO LEVEL AND SLOPE
ARE FIT. 

INPUTS:
     sigmalimit: excludes points with residuals gt sigmalimit*sigma

     xdata: the x-values at which the data points exist.
     tdata: the data points.

     zro0: the estimated constant zero offset of the data points.
     slope0: the estimated slope. 
     alpha0: the estimated coma of the first Gaussian.
     hgt0: the array of N estimated heights of the Gaussians.
     cen0: the array of N estimated centers of the Gaussians.
     wid0: the array of N estimated HPBWs of the Gaussians.

OUTPUTS:
     tfit: the fitted values of the data at the points in xdata.
     sigma: the rms of the residuals.
     zro1: the fitted zero offset.
     slope1: the fitted value of the slope.
     alpha1: the fitted coma of the first Gaussian. 
     hgt1: the array of N fitted heights. 
     cen1: the array of N fitted centers.
     wid1: the array of N fitted HALF-POWER WIDTHS.
     sigzro1: the 'error of the mean' of the fitted zero offset.
     sigalpha1: the array of errors of the sigma (coma).
     sighgt1: the array of errors of the N fitted heights.
     sigcen1: the array of errors of the N fitted centers.
     sigwid1: the array of errors of the N fitted widths.
     problem: 0, OK; -1, excessive width; -2, >50 loops; -3, negative sigmas,
    4, bad derived values.
     cov: the normalized covariance matrix of the fitted coefficients.


RELATED PROCEDURES:
   GCURV_ALLCAL (COMPUTES THE CURVE GIVEN THE PARAMETERS)
HISTORY:

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/gfit_allcal.pro)


INTENSITYCAL_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
 intensitycal_newcal
CALIBRATES INTENSITY SCALE. DOES ON SOURCE AND OFF SOURCE SCANS SEPARATELY.
INPUTS ARE:
   TCALXX, TCALYY: THE CAL VALUES FOR XX, YY.

   STOKES, THE UNCALIBRATED STOKES PARAMETERS.

and the following come from structure SCNDATA:

   SCANS, THE SET OF INDICES THAT SHOULD BE CALIBRATED. 
       FOR EXAMPLE, IN A SET OF 18 OBSERVATIONS OF CROSS PATTERN,
       YOU WANNA TO ALL OF THEM, SO SCANS=INDGEN(18).
   INDXCALON, THE INDX NR OF THE SCAN ARRAY THAT HAS THE CAL ON.
       FOR EXAMPLE, IN THE CROSS PATTERN WE USE ONLY OFF-SOURCE
       SCANS TO DO THE CALIBRATION, SO INDXCALON=[0,15]
   INDXCALOFF, THE INDX NR OF THE SCAN ARRAY THAT HAS THE CAL OFF.
       FOR EXAMPLE, IN THE CROSS PATTERN WE USE ONLY OFF-SOURCE
       SCANS TO DO THE CALIBRATION, SO INDXCALOFF=[1,14]

 polBad: int  =0 --> polA bad, 1--> polB bad else both pols ok
                makes Q 0, u,v junk
 missingcal: int   if true then there is no cal for this receiver
              set caldeflection to tsys.

OUTPUTS:
   STOKESC1, THE INTENSITY-CALIBRATED STOKES PARAMETERS.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/intensitycal_newcal.pro)


MAINBEAM_EVAL_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
 mainbeam_eval_newcal
PURPOSE: evaluate the main beam at any set of az, za (fully vectorized).
Response at beam center is UNITY TIMES THE SOURCE ANTENNA TEMPERATURE.

CALLING SEQUENCE:

   MAINBEAM_EVAL_NEWCAL, azarray, zaarray, b2dfit, mainbeam, $
   noalphaway=noalphaway

INPUTS:

   AZARRAY, ZAARRAY, the array of az and za offset from center at
which to evaluatte the beam in units of arcmin. 

   B2DFIT, the array of 2d fit coefficients described and derived in
BEAM2D_DESCRIBE. 

KEYWORDS:
   NOALPHAWAY: does NOT use the 'alpha' formulation; instead 
it uses the fourier rep of HPBW. WE DO **NOT** SET NOALPHAWAY, which
means that we do it the OLD way below.

OUTPUT: 

   MAINBEAM, the response of the main beam at the azarray, zaarray
points. 

RESTRICTIONS:

IMPORTANT NOTE: In the last statement of this pgm, notice
   that we limit the coma correction to 0.75. We do this to 
   keep the beam shape well-behaved. This is important only for
   ridiculously large values of alpha (this works up to alpha=0.5)

HISTORY:
   extracted from G2DCURV_ALLCAL.pro
   17 OCT 2000: changed definition of alpha_coma...

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mainbeam_eval_newcal.pro)


MAKE_AZZA_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
make_azza_newcal
PURPOSE: make N x N arrays of az, za values that duplicate the
 N original positions sampled in the pattern observations. 

CALLING SEQUENCE:
MAKE_AZZA, ptsperstrip, b2dfit, pixelsize, azarray, zaarray

INPUTS:

   PTSPERSTRIP, the nr of points per strip (e.g. 60)

   B2DFIT, the array of 2d fit coefficients generated by
BEAM2D_DESCRIBE (whose documentation is therein).

OUTPUTS:
   PIXELSIZE, the size of the pixels in arcmin
   AZARRAY, ZAARRAY, two arrays of azoffset, zaoffset in units of arcmin.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/make_azza_newcal.pro)


MMBMTOSTR - MOVE MUELLER BEAM PROCESSED ARRAYS TO A STRUCTURE

[Previous Routine] [Next Routine] [List of Routines]
NAME:
mmbmtostr - move  mueller beam processed arrays to a structure
SYNTAX: mmbmtostr,nrc,beamin_arr,beamout_arr,hb_arr,mm_arr
ARGS:
	nrc			  : int	index into beamin,beamout to use
	beamin_arr[M]  : holds info input to beam fits
	beamout_arr[M] : holds info output from beam fits
	hb_arr[4]      : {hdr} first header of each strip of pattern
RETURNS:
	mm			  :{} struct holding the combined info.

DESCRIPTION:
	mmtostr takes the data from the mueller matrix fitting and loads
it into a structure that is used in the archiving.


30jun01.. fixup eta mainbeam, eta sidelobe
21oct02. carl modifed to work with newcal, eliminating the intermediate
         hdr params.
09dec04.. phil update to merge this format with the old format..

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mmbmtostr.pro)


MMPROC - PROCESS 'BEETLE' SCANS (AKA HEILES SCANS)

[Previous Routine] [Next Routine] [List of Routines]
NAME:
mmproc - process 'beetle' scans (aka heiles scans)

SYNTAX:ntot=mmproc(corfile,scndata,$
          mmInfo_arr,hb_arr,beamin_arr,beamout_arr,mmCmp_arr,mmCmpChn_arr,$
          sl=sl,retsl=retsl,dirsave=dirsave,dirplot=dirplot,polBad=polBad,$
	       missingcal=missingcal,$
          tcalxx_board=tcalxx_board,tcalyy_board=tcalyy_board,han=han,$    
          board=board,rcvnam=rcvnam,sourcename=sourcename,npatterns=npatterns,$
          m_rcvrcorr=m_rcvrcorr,m_skycorr=m_skycorr,m_astro=m_astro,$ 
          mm_pro_user=mm_pro_user,$
          phaseplot=phaseplot,plot1d=plot1d,plot2d=plot2d,$
          print1d=print1d,print2d=print2d,srcprint=srcprint,
          samemm0=savemm0,onlymm4=onlymm4,$
          mm4=mm4, plt1yes=plt1yes,ps1yes=ps1yes,byChnl=byChnl,savemm4=savemm4,
          check=check,negate_q=negate_q,m7=m7,plt0yes=plt0yes,keywait=keywait,
          cumcorr=cumcorr,badScanRange=badScanRange)
 ------------------------------------------------------------------------=-
ARGS:
   corfile: name of correlator file to process path and filename.
   scndata: {scndata} information that parameterizes how scan was done. 
                      This is set by @mm0init
 ------------------------------------------------------------------------=-
RETURNS:
             ntot:long The number of beam fits done. For each pattern used
                       there will be a beamfit for each board.
 mmInfo_arr[ntot]:{mueller} hold src info, fit info for each fit.
   hb_arr[4,ntot]:{hdr} for each fit done, the hdr from the first point of 
                        each of the 4 strips.
 beamin_arr[ntot]:{beaminput} input for 2d beamfit (output of 1s fits)
beamout_arr[ntot]:{mmoutput}  output of 2d beamfit.
mmCmp_arr[ntot]:{muellerparams_carl}  if mm4 keyword set, then the
                        computed mueller matrices are returned here.
mmCmpChn_arr[128,ntot]:{muellerparams_carl}  if mm4 keyword set and /byChnl
                        keyword set then this has the mueller matrices by
                        frequency channel.
          retsl[]:{sl}  Return the scanlist from scanning the file.
                        If you reprocess the file, pass this in with sl=retsl
DEFINITIONS:
   The structures are defined in the file:
  aodefdir()/spider/pro/mm0/hdrMueller.h
It is sourced via the @mm0init call.
The  file also loads the structure {scndata} with the parameterization
of the pattern.
 ------------------------------------------------------------------------=-
KEYWORDS:
   ------------------------
   GENERAL
   ------------------------
   tcalxx_board[4]:float if negative then lookup cal values. if positive
                   then use these values for polA.brd 0 to 3.
   tcalyy_board[4]:float if negative then lookup cal values. if positive
                   then use these values for polB
     sl[]:{sl}     scanlist to use. Pass this in if you have already scanned 
                   the file.
  dirsave:string   Place the save files in this directory. The default is 
                   the current directory you are running from
  dirplot:string   Place any plotfiles in this directory. The default is 
                   the current directory you are running from
  polBad:int       if 0 then polA bad, if 1 then polB. If bad then pol will
                   be ignored (other pol copied to this one).
                   This means that Q is zero and u,v are 
                    garbage. Use when a pol is bad or not wanted. 
  missingcal :int  if set then this receiver has no cal, just use
                   avg of calon,caloff to normalize to units of Tsys.
   ------------------------
   LIMIT THE DATA TO PROCESS: default is all the scans
   ------------------------
    board: int     0..3 only process this board. The default is toprocess all
                   of the boards.
   rcvnam: string only process this receiver. The default is to process
                   all receivers. The rcvr names are:
                327,430,610,lbw,lbn,sbw,sbh,cb,xb,lbwifw,lbnifw,430ch,chlb,
                and sb750
sourcname[]: string  array of source names to process. If not provided then
                do all of the sources.
npatterns: long   if provided, then only process this many patterns
                     before returning. The default is to process all that 
                     match the above criteria.
badScanRange[2]:  long  if specified and both values > 0 then exclude
                     all scans within the range (inclusive)
excludeBrds[n] :  long exclude boards. use if bad synth or cor brd
                       brds are 0..3. -1 in first entry will ignore keyword
   ------------------------
   MUELLER CORRECTION KEYWORDS:
   ------------------------
 m_rcvrcorr:       if set, correct for rcvr,iflo. Set this if you want to
                   get mm-corrected data. If you want to check to see how well 
                   a mm matrix is applied, set this and don't set m_skycorr
                   and m_astro; then solve for the mm parametrs. they should 
                   be zero.
  m_skycorr:       if set, corrects for sky parallactic rotation. Normally
                   you set this when correcting data, and don't set it when
                   computing new matrices.
    m_astro:       if set, corrects electronics to astronomical position
                   angle and stokes V definition. The astronomical angle
                   rotation is not measured for all feeds.
mm_pro_user:string Name of alternate routine to provide mueller matrix
                   parameters. eg for lbw:mmp_lbw_08sep00_nocalcorr

   ------------------------
   PLOTTING/PRINTING KEYWORDS:
   ------------------------
  phaseplot:       if set then make phase vs freq plots for cal and source.
     plot1d:       if set then plot 1d strip fits on the screen.
     plot2d:       if set then plot 2d beam image fits on the screen.
    print1d:       if set then print 1d  fits on screen.
    print2d:       if set then print 2d  fits on screen.
   ------------------------
   WHEN COMPUTING THE MUELLER MATRIX:
   ------------------------
        mm4:       If set then compute the mueller matrix.
                   m_rcvrcorr,m_skycorr, and m_astro should not be set.
                   You need a polarized source with enough patterns to do a fit
                   in parallactic angle.
    plt1yes:       if set then plot the PA dependencies, parameters, and
                   mueller matrix elements.
     ps1yes:       if set then generate post script files of the totalpower
                   mueller matrix (same as screen output from plt1yes). An
                   example output filename would be:
                   lbw_1175_bd0_B0518+165_m4_10-AUG-2001.ps 
                   If /byChnl is set then a 2nd plot of parameter vs freq
                   will be written to:
                   lbw_1175_bd0_B0518+165_m4frq_10-AUG-2001.ps 
                   Both of these will be in the directory set by the
                   keyword dirplot=
      byChnl:      if set, compute the mueller matrix for the continuum 
                   and then for each frequency channel.
   ------------------------
   SAVING DATA IN idl SAVE FILES
   ------------------------
    savemm0:       if set then save the 1d,2d fits (hb_arr,mmInfo_arr,
                   beamin_arr,beamout_arr, and mm0procByChnl in an idl save
                   file: mm0_ddmonyy.projid.n.sav  .
                   If byChnl keyword was also set, then a second save
                   file: mm0_ddmonyy.projid.n.sav.byChnl will hold 
                   stkOffsets_chnl_arr. This is the stokes data by channel.
                   The variable mm0procByChnl=1 if the channel data was
                   also saved. These save files can be used later with the
                   /onlymm4 keyword to just do the /mm4 processing.
    savemm4:       if set and /mm4 was set then save the mueller matrix
                   computations: (mmCmp_arr,mmCmpChn_arr,mmInfo_arr,
                   mm4procByChnl) will be saved in the file:
                   mm4_ddmonyy.projid.n.sav. mmInfo_arr is a copy of what
                   would have been saved in /savemm0.
   ------------------------
   KEYWORDS FOR DEBUGGING .. normally not used
   ------------------------
    plt0yes:       if set then mm4 plots intermediate results (PA dependencies 
                   before the fit on the screen); usually not set unless
                   there are problems with the fit.
       m7:         if set, then mm4  gets the continuum using the 'M7' method.
                   It tries to excise interference from each spectrum using
                   cumcorr (it is time consuming).
    check:         checks the calc by plotting on the screen the
                   mm-corrected input data; derived MM elements should be zero.
                   Normally don't bother with this; useful for pgrm
                   development and looking into problems with fits.

 negate_q:         multiplies uncorrected xmy by -1. always use 0 here.
  keywait:         when doing the plots on the screen, wait for a keypress
                   before continuing 
 ----------------------------------------------------------------------------

DESCRIPTION:

   mmproc is the general interface to carl heiles polarization calibration
using the 'beetle' scans (rumor has it that beetles have 6 legs..). By
default it will process all receivers, sources, and boards that it finds
in the file. You can limit what is processed with the board, rcvnam, source,
etc keywords described above.

There are two basic modes for using this routine:

1. You've taken some data on a source using the beetle scans and you want
   to get the gain, Tsys, beamwidth, etc... and you are not interested in
   the stokes values: You should run this routine as:
   idl
   @mm0init
   corfile='/share/olcor/corfile.17aug02.x102.1'
   ntot=mmproc(corfile,scndata,mmInfo_arr,hb_arr,beamin_arr,beamout_arr,$
   with the optional parameters:  /phaseplot,/plot1d,/plot2,/savemm0)
   All the info you want will be returned in mmInfo_arr.
   If you are interested in the polarization info you should include:
        /m_rcvrcorr,/m_skycorr,/m_astro..   

2. You've tracked a polarized source across a large fraction of its 
   parallactic angle and you want to recompute the mueller matrix for this
   receiver. Call mmproc with:

   corfile='/share/olcor/corfile.18aug02.x102.1'
   ntot=mmproc(corfile,scndata,mmInfo_arr,hb_arr,beamin_arr,beamout_arr,
               mmCmp_arr,mmCmpChn_arr,/mm4... 
   with optional keywords: /phaseplot,/plot1d,/plot2,/savemm0, 
               /plt1yes,/ps1yes,/savemm4
   The mueller matrix info will be in mmCmpm_arr.

   In either case, use dirplot,dirsave to set the directory for the
   output save files or ps files (if they were requested).

After running the routine, you can replot the stripfits or 2d beam images
with the routines: plot_beam1d, or plot_beam2d

The startup of idl is:
   idl
   @mm0init
   corfile='/share/olcor/...'
   ntot=mmproc(corfile...)
   You can then continue calling mmproc as many times as you want.

SEE ALSO:
   mmplot_beam1d,mmplot_beam2d,mm_structureDefs

********************BEGIN WARNING NR 2*****************************

   IF YOU RUN THIS WITHOUT MUELLER CORRECTING THE INPUT DATA,
THEN THE BEAM MAPS FOR THE POLARIZED STOKES PARAMETERS WILL NOT BE FOR
THE ***REAL*** STROKES PARAMETERS, BECAUSE WILL NOT HAVE BEEN CALIBRATED!
*********************END WARNING NR 2*****************************

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mmproc.pro)


MMPROCLIST - DO MUELLER PROCESSING FOR A SET OF DATA FILES

[Previous Routine] [Next Routine] [List of Routines]
NAME:
mmproclist - do mueller processing for a set of data files

SYNTAX: ntot=mmproclist(listfile,scndata,$
           mmInfo_arr,hb_arr,beamin_arr,beamout_arr,mmCmp_arr,mmCmpChn_arr,$
          fnameind=fnameind,maxfiles=maxfiles,_extra=_e
ARGS:
   listfile: string. filename holding names of correlator files 
                     (one per line). semi-colon as the first char is a comment.
   scndata: {scndata} information that parameterizes how scan was done.
                      This is set by @mm0init
KEYWORDS:
 fnameind: int       index in line for filename (count from 0)
 maxfiles: int       maximum number of files to process (default is
					  99999L)
   _extra: {}      parameters that can be used by mmproc. Any parameters
                   will be applied to all of the files.
RETURNS:
             ntot:long The number of beam fits done. For each pattern used
                       there will be a beamfit for each board.
 mmInfo_arr[ntot]:{mueller} hold src info, fit info for each fit.
   hb_arr[4,ntot]:{hdr} for each fit done, the hdr from the first point of
                        each of the 4 strips.
 beamin_arr[ntot]:{beaminput} input for 2d beamfit (output of 1s fits)
beamout_arr[ntot]:{mmoutput}  output of 2d beamfit.
mmCmp_arr[ntot]:{muellerparams_carl}  if mm4 keyword set, then the
                        computed mueller matrices are returned here.
mmCmpChn_arr[128,ntot]:{muellerparams_carl}  if mm4 keyword set and /byChnl
                        keyword set then this has the mueller matrices by
                        frequency channel.
          retsl[]:{sl}  Return the scanlist from scanning the file.
                        If you reprocess the file, pass this in with sl=retsl
DESCRIPTION:
   Call mmproc for every filename in listfile. Return all of the 
data in one array of structures. An example of the listfile is:

;  a comment
/share/olcor/calfile.19mar01.a1400.1
/share/olcor/calfile.20apr01.a1446.1
/share/olcor/calfile.20apr01.a1489.1
/share/olcor/calfile.20mar01.a1389.1

It would process all 4 files and return the data in mm,mmInfo.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mmproclist.pro)


MMPROC_SETUP - SETUP ROUTINE FOR MMPROC

[Previous Routine] [Next Routine] [List of Routines]
NAME:
mmproc_setup - setup routine for mmproc
SYNTAX: numPat=mmproc_setup(corfile,scndata,lun,retsl,indForPat,$
                   brdToUse,nfits,
                  mmInfo_arr,hb_arr,beamin_arr,beamout_arr,stkOffsets_chnl_arr,
                   filesavenamemm0, filesavenamemm4,onlymm4=onlymm4,$
                  tcalxx_board=tcalxx_board,tcalyy_board=tcalyy_board,$
                   sl=sl,board=board,rcvnam=rcvnam,byChnl=byChnl,
                   sourcename=sourcename,npatterns=npatterns,
                   dirsave=dirsave,dirplot=dirplot,badScanRange=badScanRange,$
                   excludeBrds=excludeBrds)
ARGS:
   corfile: string filename with correlator data (includes path)
   scndata: {scndata} observing pattern parameters.
RETURNS:
           numPat: int total number of patterns to do.
          retsl[]:{sl} scan list of entire file
 indForPat[numPat]: int pointers into retsl[] for start of each pattern to do.
          brdToUse: int   -1--> all boards, or board number 0..3
             nfits: int    total number of fits to do. NumPat*numbrds
  mmInfo_arr[nfits]:{mueller} hold results from fits
 beamin_arr[nfits]:{beaminput} input data before fits
beamout_arr[nfits]:{mmoutput} output from fits
   filesavenamemm0: string save filename for mm0 1d fits
   filesavenamemm4: string save filename for mm4 mueller matrix computations

KEYWORDS:
tcalxx_board[]: float if positive then use these value for calxx 
tcalyy_board[]: float if positive then use these values for calyy
   sl[]:{sl}   scanlist user passes in (optional)
  board: int   0..3 process only this board
 rcvnam: string process only this receiver
 sourcename[]:string process only these sources
 npatterns:int stop after processing this many patterns.
   dirsave:string  directory for save files. default is current
   dirplot:string  directory for plot files. default is current
   onlymm4:    if set then only do mm4 processing. input the save file
               and do the subsetting they want
 badScanRange[2]; long  if present and values >0 then exclude any scans between these to limits(inclusive).

DESCRIPTION: 
   Scan the file, figure out where the patterns are, make any subset
of file they request, initializes variables.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mmproc_setup.pro)


MM_CORR

[Previous Routine] [Next Routine] [List of Routines]
NAME:
mm_corr
PURPOSE: Mueller-correct the stokes data .stkoffsets_chnlxx 
   (data that have already been calibrated with the cal). 

INPUTS:

   MM_PRO_USER: String; If specfied, use this proc to determine the
mm coefficients instead of the default.  Set equal to '' for default. 

   cfr: float  center frequency in Mhz for the data.
   stokesc1[128,4,242] phase and intensity calibarated data to correct
   
   azencoders[60,4]: az encoder for the strip points
   zaencoders[60,4]: az encoder for the strip points
           stkoffsets_chnl[128,4,60,4,*] and stkoffsets_chnl_cal[128,4,2,*].

KEYWORDS:

   M_RCVRCORR is normally set, corrects for rcvr--this is the
matrix that we put so much effort into calibrating!

   M_SKYCORR is normally set, corrects for sky rotn wrt az arm 
during tracking.

   M_ASTRO is normally set, rotates PA's to astronomical definition.
this rotation is not known for all rcvrs. at the moment it is known only
for LBW and, maybe, CB. It is defined in the mmp.*** files.

OUTPUTS:
   BEAMIN: structure that contains the mm-corrected stkoffsets_chnlxxx data
(the correction is done in place)

   MCORR: the first three digits of this integer number tell
whether: first digit, m_rcvr was applied; second, m_ksy was applied;
third, m_astro was applied.

NOTE:
   the data stokesc1 will be corrected according to the keywords
   m_rcvrcorr,m_skycorr,m_astro
   the cal data stokesc1[*,*,0:1] will  only be corrected for
   m_rcvrcorr no matter what the other keywords are set to. This will
   let you see how well the rcvr_cor is doing since all of the polarized
   power should end up in U for the cal difference.
   Because of this we can use beamin_arr.azenc.. since we only need
   the 240 locations for the paralactic angle correction, not the cals.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mm_corr.pro)


MM_CORR_STOKESC1

[Previous Routine] [Next Routine] [List of Routines]
NAME:
mm_corr_stokesc1
PURPOSE: mueller-correct the stokesc1 array IN PLACE. This is the primary
data generated by CROSS3_GENCAL.

CALLING SEQUENCE:

   MM_CORR_STOKESC1, m_tot, m_astron, azmidpntn, zamidpntn,$
   m_rcvrcorr= m_rcvrcorr, m_skycorr= m_skycorr, m_astro=m_astro


INPUTS:

       M_TOT, the mueller matrix for the electronics system as defined
in the AOTM

       M_ASTRON, the mueller matrix to convert to astronomy as defined
in the AOTM

   stokesc1[128,4,242] - intenisty and phase calibrated data
   azencoders[60,4]    - az pos for non-cal data
   zaencoders[60,4]    - za pos for non-cal data


KEYWORD:

   M_RCVRCORR: setting it includes the M_TOT correction.

   M_SKYCORR, setting it includes the sky correction; default is zero.

   M_ASTRO: setting includes the M_ASTRON correction

OUTPUTS:

   stokesc1: the corrected data.

NOTE:
   the stokesc1 data will be corrected according to the keywords
   m_rcvrcorr,m_skycorr,m_astro
   the data stokesc1[128,4,0:1] (the cals) will only be corrected for
   m_rcvrcorr no matter what the other keywords are set to. This will
   let you see how well the rcvr_cor is doing since all of the polarized
   power should end up in U for the cal difference.
   az,za encoders are 60,4 and are the positions for the strip data
   (not including the cals)

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/mm_corr_stokesc1.pro)


PA1FIT

[Previous Routine] [Next Routine] [List of Routines]
NAME:
pa1fit
 PURPOSE: LEAST-SQUARE FITS 

   INT = A*COS(PA1) + B*SIN(PA1) + C (C is optional)

DISCARDS POINTS EXCEEDING 3 SIGMA!!!

 CALLING SEQUENCE:

PA1FIT, pa1, const, int, $
   coeffs, coeffsp, sigma, resids, niterations, cov

 INPUTS:

   PA1, the array of angles. Units are DEGREEES
   CONST: set equal to zero to NOT fit the C above, nonzero to fit
the C. 
   INT, the array of intensities.

KEYWORD:
   NODISCARD: setting it prevents testing residuals and discarding
bad points.

 OUTPUTS:
   COEFFS = fltarr( 3,2):
       COEFFS[*,0] are the fitted coefficients
       COEFFS[*,1] are the errors

   COEFFSP = fltarr( 2,2). 
       COEFFSP[0,*] is the amplitude and its error.
       COEFFSP[1,*] is the phase angle and its error, units 
       DEGREES!!!
   
   SIGMA: the sigma of the residuals from the fit
   RESIDS: the array of residuals
   NITERATIONS: number of iterations used in discarding 3sig
points.
   COV: the normalized covariance matrix

 EXAMPLE:

 MODIFICATION HISTORY:

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/pa1fit.pro)


PHASECAL_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
 phasecal_newcal
PURPOSE: Least squares fit for the phase of the correlated cal versus
frequency, and correct the data if requested by CORRECTOPTION.

INPUTS: 
   FRQ[ 2048]: the aray of frequencies, centered at the frequency
at which you want the zero of phase to be calculated. Units are MHz.

   DPDF (from SCNDATA_DPDF): the initial guess for (dphase/dfreq),
radians/MHz. (from structure SCNDATA)

   STOKESC1[2048, 4, *]: the input Stokes parameters

   PHASECHNLS: the particular chnls to include in the fit (from 
structure SCNDATA).

   INDXS: the set of indices for the third elements of STOKESC1 to
correct, if correction is chosen as the option.

   INDXCALOFF: the set of indices in the third element of STOKESC1
to use as caloff for the calibration.

   INDXCALON: the set of indices in the third element of STOKESC1
to use as calON for the calibration.

KEYWORD: CUMCORR excises interference from the calon/caloff spectra

   CORRECTOPTION: 
       if 1, it corrects STOKESC1; 
       if 0, it does no correction;
       if -1, it corrects for the slope but not the phase.


OUTPUT:

   if CORRECTOPTION is 1, then STOKESC1 are returned as
phase-calibrated stokes parameters.

   PHASE_OBSERVED is the array of observed phases of the cal
deflection. 

   OZERO and OSLOPE are the linear fit coefficients and errors,
the units are **** RADIANS ***** and ***** RADIANS/MHZ *****

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/phasecal_newcal.pro)


PHASEFIT_F

[Previous Routine] [Next Routine] [List of Routines]
NAME:
   PHASEFIT_F

PURPOSE:
    Fit real and imag simultaneously to the function
   CA*cos(KK*frq) - SA*sin(KK*frq) +
   SA*cos(KK*frq) + CA*sin(KK*frq) = [realpart, imagpart]
    and solves for CA, SA, and KK.

    A three step process:
   1. A linear fit for aareal, etc, using the guess for dpdf
   2. A nonlinear fit for dpdf using the above values for aareal, etc.
   3. A simultaneous onlinear fit for all five parameters.

    In each step, bad points are tested for and discarded using the
   criterion residual gt sigma3 * sigma.

CALLING SEQUENCE:
    PHASEFIT_A, frq, realpart, imagpart, dpdf, coeffs, sigcoeffs, sigma

INPUTS:
     frq: freq, in MHz, of the XY products; offset should be
   subtracted so that the values cluster around zero.
     realpart: the real part of the XY products. 
     imagpart: the imaginary part of the XY products. 
     dpdf: the ESTIMATE for d(phase)d(freq), i.e. the estimate for KK.
       THE UNITS FOR DPDF ARE DEGREES PER MHZ.
       THE PROGRAM SOLVES FOR KK, WHICH IS IN RADIANS PER MHZ.

OUTPUTS:
     coeffs: = [CA, SA, KK*!radeg]. 
       NOTE: COEFFS[2] IS IN DEGREES/MHZ, WHILE KK IS IN RAD/MHZ.
     sigcoeffs: the uncertainties of the coefficients.
     sigma: the mean error of the residuals.
     TFIT: the fitted datapoints.
     APHI is the phase angle at frq=0, units degrees. Thus:
   the DATA have PHASEANGLE = ATAN(IMAGPART,REALPART) radians
       or  !radeg * ATAN(IMAGPART,REALPART) degrees
   the FIT has   PHASEANGLE = APHI + COEFFS[2]*FRQ degrees
     SIGAPHI is the error in APHI
     PROBLEM is nonzero if there is a problem...value depends on the problem

RESTRICTIONS:
    As with any nonlinear least squares fit, the estimate should
   be pretty good. Zero estimates for dpdf are replaced intgernally
   by dpdf=20 degrees over the full bandpass.


RELATED PROCEDURES:
   PHASEGCURV

HISTORY:
   Written by Carl Heiles. 19 FEB 1999.
   Corrected and updated 4 MAR 2000.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/phasefit_mar01.pro)


PLOT_BEAM1D

[Previous Routine] [Next Routine] [List of Routines]
NAME:
plot_beam1d
PLOT THE RESULT OF THE 1-D STRIPFITS. 4 PANELS, ONE FOR EACH STOKES
PARAMETER; 4 COLORS, ONE FOR EACH STRIP.

INPUTS:
   B_0, the structure containing frequency info and scan nr. This
is used only to annotate the plot title.

   BOARD, the board number. Again, used only to annotate the plot
title. 

   TOTOFFSET, the total angular offsets of the strip positions in
arcmin. totoffset[60,4]: 60 point per strip, 4 strips. units are ARCMIN

   STOKESOFFSET_CONT, the stokes parametes for the strips in
Kelvins. stokesoffset_cont[ 4, 60, 4]: 4 stokes, 60 pts/strip, 4 strips.

   TEMPFITS, the ls fits to the data in stokesoffset_cont.
tempfits[ 4, 60, 4], as in stokesoffset_cont.

COMMON:
   uses the colors defined in Carl's PLOTCOLORS common block.

OUTPUTS: THE PLOT IN WINDOW 31. THE WINDOW IS CREATED IF IT IS NOT
THERE.

NOTE: this will do any nr of strips up to four; it selects the number
according to the third dimensional element of tempfits.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/plot_beam1d.pro)


PLOT_BEAM2D

[Previous Routine] [Next Routine] [List of Routines]
NAME:
plot_beam2d
PURPOSE: display the 2d beam; return the total beam map optional

INPUTS:

   BEAMOUT, the structure containing the output data.

   WINDOWSIZE, the size of windows to give the image in pixels. The
windows are square and if the existing windoes are not of the proper
size they are regenerated. We use windows 0 and 1.

OUTPUTS:

OPTIONAL OUTPUTS:

   TOTALBEAM_MAG, the windowsize by windowsize image of the total
beam that is normally displayed on the little 200 by 200 window in the
bottom right corner.

KEYWORDS:

   SHOW: if set, it gives images of the fitted beam and sidelobe

   NTERMS, the number of terms to use in reconstructing the
sidelobe structure from the Fourier fit to the sidelobe structure. Use
of NTERMS=8 forces the reconstruction to duplicate the measured points
because there were 8 measured data points (2 on each end of the strip,
so they are located every 45 deg around the circle). Use of NTERMS<8 is
equivalent to least squares fitting those 8 points with a Fourier series
having the chosen value of NTERMS, which means that the fitted curve
will not go through the points. Seeing as the 'measured points' are really
the result of model fits to the sidelobes with Gaussians, and that these
fitted parameters are not perfect, I don't believe using NTERMS=8 is
justified. Rather, some smaller value is better. Currently I am using 6.
This needs some investigation by looking at the sidelobe coefficients,
the original strip data, and the reconstructed results using different
values of NTERMS to see what's appropriate. 

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/plot_beam2d.pro)


POLBEAM_EVAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
polbeam_eval
PURPOSE: evaluate the POLARIZED STOKES PARAMETER MAIN BEAM at any set
of az, za (fully vectorized). 

INPUTS:

   ARCMIN, the hpbw specified in the original pattern observation,
arcmin. 
   AZARRAY, ZAARRAY, the array of az and za in units of ARCMIN.
   B2DFIT, the array of 2d fit coefficients described and derived in
BEAM2D_DESCRIBE. 

   NSTK, the Stokes parameter number: 1,2,3 for Q,U,V.  Don't set
NSTK=0

OUTPUTS: 
   SQUINTBEAM, the squint portion of the polarized beam evaluated
at (azarray, zaarray)
   SQUASHBEAM, the squint portion of the polarized beam evaluated
at (azarray, zaarray)

RESTRICTIONS: Some IMPORTANT NOTES: see documentation for
MAINBEAM_EVAL.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/polbeam_eval.pro)


PRINT_BEAM1D - PRINTS OUT THE 1-D BEAMFITS.

[Previous Routine] [Next Routine] [List of Routines]
NAME:
print_beam1d - PRINTS OUT THE 1-D BEAMFITS.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/print_beam1d.pro)


SIDELOBE_EVAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
sidelobe_eval
PURPOSE: evaluate the sidelobes at an arbitrary set of az, za given the
Fourier coefficients describing hgt, cen, and wid.

CALLING SEQUENCE:
SIDELOBE_EVAL, nterms, fhgt, fcen, fhpbw, az, za, sidelobe

INPUTS:

   NTERMS, the number of terms to use in evaluating the inverse
FFT. For our application with 4 strips, we use NTERMS=6).

   FHGT, FCEN, FWID, the set of complex Fourier coefficients for
the height, center, and width of the Gaussians that describe the
sidelobes.

   AZ,ZA: the position at which to evaluate the sidelobes. Units
are ARCMIN. 

OUTPUT:    

   SIDELOBE, the amplitude of the sidelobe.

NOTE ON UNITS FOR THE WIDTH: The input width is in in units of HPBW. In
the evaluation, we use 1/e widths, which is why we multiply by 0.6...

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/sidelobe_eval.pro)


SQ2DFIT_ALLCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
sq2dfit_allcal
   SQFIT_ALLCAL, nrstrip, nstk, xdata, stokesoffset_cont, $
        tfit, sigma, stripfit, sigstripfit, problem, cov

PURPOSE:
FIT    (1) INTENSITIES OF CENTRAL BEAM AND THE TWO SIDELOBES,
       (2) BEAM SQUINT (FIRST DERIVATIVE OF GAUSSIAN) AND
       (3) SQUASH (DIFF IN BEAMWIDTHS) FOR ALL POL STOKES AND ONE STRIP.
   note: added squash_avg to ls soln 31jul2003.
       For deriving polarization leakage, beam squint and
       beam squash. the parameters
       of the original gaussian are already known.

CALLING SEQUENCE:
    SQFIT_ALLCAL, nrstrip, nstk, offset, stokesoffset_cont, $
        tfit, sigma, stripfit, sigstripfit, problem, cov

INPUTS:
     xdata: the x-values at which the data points exist.
     STOKESOFFSET_CONT, the datapoints
     STRIPFIT[ *, 0, NRSTRIP]: THE XPY POWER FOR THE STRIP (ALREADY DERIVED)

OUTPUTS:
     tfit: the fitted values of the data at the points in xdata.
     sigma: the rms of the residuals.
     STRIPFIT[ *, NSTK, NRSTRIP] coeffs: array of fitted parameters for
       this stokes param NSTK, this strip NRSTRIP
     SIGSTRIPFIT: the array of errors of coeffs.
     problem: 0, OK;  -3, negative sigmas; 4, bad derived values.
     cov: the normalized covariance matrix of the fitted coefficients.

INTERMEDIATE RESULTS...
       COEFFS: the array of fitted parameters, which are:
       [zero offset, main beam polarized intensity, squint, squash,
               left sidelobe polarized intensity,
               right sidelobe polarized intensity]

RELATED PROCEDURES:
       GCURV_ALLCAL
HISTORY:
       Derived from sqfit.pro, 13sep00
   ch eliminated subtracting the polarized tsys offsets. 04jan03
       2003feb06: factor of two error in squint, squash fixed. earlier
squint/squash were factor of two too small. In PASP article, see
equation 20: the factor (I/2) sits iin front of the second term instead
of the factor (I). this factor of two repair is put in at the end of
this program.
   2003jun31: included squoosh in the fit, so the squash 
is now like a fourier series with 0 and 2pa terms. the output array
   2003aug10: fixed squash, which has been faulty forever, and squoosh.
added keyword SQUOOSH

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/sq2dfit_allcal.pro)


SQFIT_NEWCAL

[Previous Routine] [Next Routine] [List of Routines]
NAME:
sqfit_newcal

   SQFIT_NEWCAL, nrstrip, nstk, offset, tdata, $    
        tfit, sigma, stripfit, sigstripfit, problem, cov

PURPOSE: Given a single strip scan, the previously accomplished Stokes
I fit to that scan with XPYFIT_NEWCAL, and a single one of the three
polarized Stokes parameters, fit:
       (1) Intensities of polarized power in the central beam and the 
two sidelobes (the "gain error")'
   (2) Beam squint (first derivative of Gaussian) and 
   (3) Squash (diff in beamwidths) for ALL POL STOKES AND ONE STRIP.
   For deriving polarization leakage, beam squint and 
   beam squash. the parameters
   of the original gaussian are already known.

CALLING SEQUENCE:
    SQFIT_NEWCAL, nrstrip, nstk, offset, tdata, $    
        tfit, sigma, stripfit, sigstripfit, problem, cov

INPUTS:
       NRSTRIP, the sequential number of the strip (0 to NRSTRIPS-1)

   NSTK, the Stokes parameter. (1,2,3) = (Q,U,V).

       OFFSET, the angular offset from the assumed center IN UNITS OF
THE ASSUMED HPBW.

   TDATA, the Stokes Q, U, or V system temperatures at each point. 
Units are Kelvins. 

   STRIPFIT, the input/output data array.  STRIPFIT[ *, 0, NRSTRIP] are
the PREVIOUSLY-DERIVED results for Stokes I

OUTPUTS:
     TFIT: the fitted values of the data at the points in offset.

     SIGMA: the rms of the residuals (TDATA- TFIT).

     STRIPFIT[ *, NSTK, NRSTRIP], SIGSTRIPFIT: array of fitted
parameters and errors for this stokes parameter NSTK, this strip
NRSTRIP. See full documentation in BEAM1DFIT.PRO .

     SIGSTRIPFIT: the array of errors of coeffs.

     problem: 0, OK;  -3, negative sigmas; 4, bad derived values.

     cov: the normalized covariance matrix of the fitted coefficients.

INTERMEDIATE RESULTS...
   COEFFS: the array of fitted parameters, which are:
   [zero offset, main beam polarized intensity, squint, squash, 
       left sidelobe polarized intensity,
       right sidelobe polarized intensity]

RELATED PROCEDURES:
   GCURV_ALLCAL
HISTORY:
   Derived from sqfit.pro, 13sep00
   NEWCAL mod by carl, 9 oct 2002
   2003feb06: factor of two error in squint, squash fixed. earlier
squint/squash were factor of two too small. In PASP article, see
equation 20: the factor (I/2) sits iin front of the second term instead
of the factor (I). this factor of two repair is put in at the end of
this program.

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/sqfit_newcal.pro)


XPYFIT_NEWCAL

[Previous Routine] [List of Routines]
NAME:
xpyfit_newcal
PURPOSE:
      For each 1-d scan (strip), and given the total intensity Stokes I, 
fits 3 Gaussians, two to sidelobes and one to main beam, for 
the strip. The main beam Gaussian has a skew parameter to represent coma.

INPUTS:
   NRSTRIP, the sequential number of the strip (0 to NRSTRIPS-1).

   OFFSET, the angular offset from the assumed center IN UNITS OF
THE ASSUMED HPBW. 

   PWR, the Stokes I total system temperature at each
observed point. Units are Kelvins.

OUTPUTS:
   TFIT, the fitted line to the datapoints.

   SIGMA, the rms of the residuals (PWR- TFIT)

   STRIPFIT, SIGSTRIPFIT, which are defined in the documentation for
BEAM1DFIT.PRO

   PROBLEM: indicates a problem with the nonlinear Gaussian fits.
I'm not sure that this is properly defined.

   COV, the covariance matrix, or its normalized counterpart, in the
nonlinear Gaussian fits. I'm not sure this is properly defined.

NOTES;
   SIGMALIMIT, is the limit beyond which points are discareded, is
hardwire defined here. 

   ALL INPUT ANGULAR WIDTHS, OFFSETS and ALL OUTPUT ANGLES are in units
of the original guessed HPBW.  They are not 1/E and the units are NOT
arcmin. 

HOW IT WORKS: 
   FIRST STAGE: It tries three Gaussians, one for each sidelobe and
one for the main beam. It checks the width and errors
for the sidelobe fits, and if they ok it returns. If one Gaussian has
parameters that are out of bounds, it goes to the second stage.

   SECOND STAGE: It tries to fit two Gaussians, one for a sidelobe and
one for the main beam.  There are two possible cases: one with the
sidelobe on the left, one with the sidelobe on the right.  It checks
both solutions and picks the one it likes best, based on the same
parameter test used in the first stage.  If neither is acceptable, it
goes to the THIRD STAGE. 

   THIRD STAGE: It fits a Gaussian to the main beam only--no sidelobes.

modified 28 oct 2000 to check widthmin as well as widthmax.
modified 19 oct 2002 for newcal

(See /pkg/rsi/local/libao/phil/spider/pro/mm0/xpyfit_newcal.pro)