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)
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)
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)
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)
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)
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)
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
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)
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)
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)
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)
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)
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)
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)
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)
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)
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
OUTPUTS:
STOKESC1, THE INTENSITY-CALIBRATED STOKES PARAMETERS.
(See /pkg/rsi/local/libao/phil/spider/pro/mm0/intensitycal_newcal.pro)
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)
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)
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)
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
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)
------------------------------------------------------------------------=-
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.
------------------------
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.
------------------------
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)
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)
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)
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
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)
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)
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)
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)
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)
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)
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)
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)
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)
NAME: print_beam1d - PRINTS OUT THE 1-D BEAMFITS.
(See /pkg/rsi/local/libao/phil/spider/pro/mm0/print_beam1d.pro)
NAME:
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)
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)
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)
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)
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)