# automap.g # Author: Steven T. Myers (NRAO) 2001-12-17 (revised) # Set of functions to process VLA continuum polarization data # # Contents: # makemsf(par) run fitstoms # automap(par) calibrate and image continuum poln data # killant(msf,ants) delete antennas using flagger # quackall(msf) quack first 10s from each scan # atan2(x,y) output atan in range -180 to +180 # regang(x) regularize degrees to lie from -180 to +180 # # Usage: # # 1. define parameters in record par # 2. if needed, make ms from fits using makemsf(par) # 3. run automap(par) # # This script is designed to work on point sources or only slightly # to moderately resolved sources (though the calibrators 0137+331 and # 1331+305 are significantly resolved, and seem to work OK). # # vlagaincor.g # Author: Steven T. Myers (NRAO) 2001-12-17 (revised) # Set of functions to apply opacity and gain correction to vla data # # Contents: # vlagaincor(par) apply gain correction to DATA column # opac(tau,el) calculate opacity correction at vector el # gaincv(a1,a2,el) gain curve for vecs an1,an2,el # vlagainfill(band) fill gain coefficient array # # Usage: # # 1. define parameters in record par # 2. if needed, make ms from fits using makemsf(par) from automap.g # otherwise BE SURE TO HAVE A BACKUP COPY OF THE MS! # 3. run vlagaincor(par) # include 'imager.g'; include 'calibrater.g'; include 'image.g'; include 'ms.g'; include 'misc.g'; include 'flagger.g'; include 'sysinfo.g'; include 'measures.g'; include 'table.g'; include 'quanta.g'; # # Some setup variables - put into record par # THESE ARE EXAMPLES - REPLACE WITH VALUES APPROPRIATE FOR YOUR DATA # par.pwdir := '/home/sandrock/smyers/Testing/2001-12-12/'; par.prefix := '20010204.kband.'; # Prefix for files generated par.fitsf := '20010204AK.UVF'; # Name of UVFITS file par.msf := '20010204AK.ms'; # Name of ms par.abscal := '1331+305'; # Name of primary flux cal par.absid := 6; # FIELD_ID of primary flux cal par.angcal := '1331+305'; # Polarization angle calibrator par.absang := 66.00; # R-L Phase Difference of PCal par.polcal := ['1159+292']; # Name of d-term calibrator par.polid :=1; # FIELD_ID of d-term calibrator par.refsrc := ['0713+438', '0854+201', '0927+390', '1159+292', '1229+020', '1256-057', '1310+323']; # Names of secondary sources par.nfld := 8; # Number of targets par.fieldlist := ['0713+438', '0854+201', '0927+390', '1159+292', '1229+020', '1256-057', '1310+323', '1331+305']; # Names of targets par.gsolint := 60.0; # solution interval (secs) for gain par.grefant := 8; # reference antenna for gain sols par.clniter := 200; # Number of clean iterations desired #par.clcell := '0.4arcsec'; # image cell size (Xband) par.clcell := '0.04arcsec'; # image cell size (Kband) par.imsize := 256; # image size in pixels # # paramters for gain correction # par.tau0 := 0.045; par.band := 'k'; # # Now define functions that use the parameters # # Make the ms from fits # makemsf := function(par) { # # Input record parameters: par.fitsf FITS file name (string) # par.msf MS file name (string) # m:=fitstoms(fitsfile=par.fitsf, msfile=par.msf); m.summary(); m.done(); print 'Data loaded into', par.msf; } # # The automapping function itself, use par record defined above # automap := function(par) { # # Input record parameters: par.pwdir Working directory name (string) # par.prefix Prefix for files generated (string) # par.msf Name of MS file (string) # par.abscal Name of primary flux cal (string) # par.absid FIELD_ID of primary flux cal (int) # par.angcal Polarization angle calibrator (string) # par.absang R-L Phase Difference of PCal (float) # par.polcal Name of d-term calibrator (string) # par.polid FIELD_ID of d-term calibrator (int) # par.refsrc Names of secondary sources (str vec) # par.nfld Number of total targets (int) # par.fieldlist Names of all targets (str vec) # par.gsolint Solution interval (s) for gain (float) # par.grefant Reference antenna for gain sols (int) # par.clniter Number of clean iterations (int) # par.clcell Image cell size (string) # par.imsize Image size in pixels (int) # msf := par.msf; # # Assume ms has been created - check with default os tool # if( dos.fileexists(file=msf) ) { print 'Found MS file ',msf; } else { fail "MS file does not exist"; } msfilename := msf; # # Copy setup variables from par into internal variables # prefix := par.prefix; # File prefix pwdir := par.pwdir; # Working directory abscal := par.abscal; # Name of primary flux cal absid := par.absid; # FIELD_ID of primary flux cal angcal := par.angcal; # Polarization angle calibrator absang := par.absang; # R-L Phase Difference of PCal polcal := par.polcal; # Name of d-term calibrator polid := par.polid; # FIELD_ID of d-term calibrator refsrc := par.refsrc; # Names of secondary sources nfld := par.nfld; # Number of targets fieldlist := par.fieldlist; # Names of targets gsolint := par.gsolint; # solution interval (secs) for gain grefant := par.grefant; # reference antenna for gain sols clniter := par.clniter; # Number of clean iterations desired clcell := par.clcell; # Cell size imsize := par.imsize; # Image size in pixels # # Basic editing # # fg:=flagger(msfilename); # fg.setflagmode('unflag'); # fg.query(query='FIELD_ID >= 0'); # fg.setflagmode('flag'); # fg.filter(column='DATA', operation='range', comparison='Amplitude', # range=['0.00001Jy','1.34Jy']); # fg.done(); # print 'Data flagged'; # # Set default source model for 3C286 # Note that the MODEL_DATA columns now default to 1.0 when imager starts # But reset them all explicity in case you are re-running this function # imagr:=imager(msfilename); imagr.setjy(); # # Solve for G Jones (amp & phase selfcal) # caltablename:=spaste(prefix, 'gcal1'); cal:=calibrater(msfilename); cal.setdata(); # First apply parallactic angle correction (necessary) cal.setapply(type='P', t=5.0); # Solve for gain, solution interval = gsolint cal.setsolve(type='G', t=gsolint, refant=grefant, table=caltablename); cal.solve(); print 'Gain calibration complete'; # # Establish the flux density scale # outtablename:=spaste(prefix, 'gcalout1'); cal.fluxscale(tablein=caltablename, tableout=outtablename, reference=abscal, transfer=refsrc); print 'Flux density calibration determined from ',abscal; # # Correct the data cal.setapply(type='P', t=5.0); cal.setapply(type='G', t=0.0, table=outtablename); cal.correct(); cal.done(); print 'Gain and Flux density calibration applied'; # # Loop over IFs # for ( spid in 1:2 ) { # # Make a map of primary calibrater, natural weighting # ifstr := sprintf('.if%d',spid); calprefix:=spaste(prefix, 'abscal', ifstr); modfilename:=spaste(calprefix, '.model'); clnfilename:=spaste(calprefix, '.restored'); resfilename:=spaste(calprefix, '.residual'); psffilename:=spaste(calprefix, '.psf'); imagr.setdata(fieldid=polid, mode='channel', nchan=1, start=1, step=1, spwid=spid); imagr.setimage(nx=imsize, ny=imsize, cellx=clcell, celly=clcell, stokes='IQUV', fieldid=polid, start=1, step=1, nchan=1, mode='channel', spwid=spid); imagr.make(modfilename); # imagr.weight('natural'); # imagr.makeimage(type='psf',image=psffilename) # bmaj:=0.; bmin:=0.; bpa:=0.; # imagr.fitpsf(psf=psffilename, bmaj=bmaj, bmin=bmin, bpa=bpa) imagr.clean(algorithm='clark', niter=clniter, gain=0.1, threshold='0.0Jy', model=modfilename, image=clnfilename, residual=resfilename); # # Extract peak (I,Q,U,V) # img1:= image(modfilename); pmax:= array(0.0, 4); for (plane in 1:4) { pregion:= drm.box([1,1,plane], [imsize,imsize,plane]); img1.statistics(statsout=pstats, region=pregion); if (abs(pstats.min) > abs(pstats.max)) { pmax[plane]:= pstats.min; } else { pmax[plane]:= pstats.max; }; }; img1.done(); # ival:= pmax[1]; qval:= pmax[2]; uval:= pmax[3]; # # Set the polarization model, applying the AIPS PCAL constraint # imagr.setjy(fieldid=polid, spwid=spid, fluxdensity=[ival, qval, uval, 0.0]); printf('IF%d Polarization model computed for %s\n',spid,polcal); # # Solve for D Jones from polarization calibrator # poltablename:=spaste(prefix, 'dcal1'); cal:= calibrater(msfilename); # Data selection parameters fldid:=polid selectr := sprintf('FIELD_ID==%d && SPECTRAL_WINDOW_ID==%d',fldid,spid) cal.setdata(msselect=selectr); # Apply parallactic angle and previous gain corrections cal.setapply(type='P', t=5.0); cal.setapply(type='G', t=0.0, table=outtablename); # Solve for single day-long D-terms, with 10min pre-avg time cal.setsolve(type='D', t=86400.0, preavg=600.0, table=poltablename); cal.solve(); # # Correct all other sources for G Jones and D Jones # selectr := sprintf('SPECTRAL_WINDOW_ID==%d',spid) cal.setdata(msselect=selectr); cal.setapply(type='P', t=5.0); cal.setapply(type='G', t=0.0, table=outtablename); cal.setapply(type='D', t=0.0, table=poltablename); cal.correct(); cal.done(); printf('IF%d Polarization model applied to all sources\n',spid); # # Make a polarization-calibrated (I,Q,U,V) image for all sources # fldtablename:=spaste(msfilename, '/FIELD'); fldtab:= table(fldtablename); fldname:= fldtab.getcol('NAME'); fldtab.close(); for (ifld in 1:nfld) { imagr.setdata(fieldid=ifld, mode='channel', nchan=1, start=1, step=1, spwid=spid); imagr.setimage(nx=imsize, ny=imsize, cellx=clcell, celly=clcell, stokes='IQUV', fieldid=ifld, start=1, step=1, nchan=1, mode='channel', spwid=spid); srcprefix:= spaste(prefix, fldname[ifld], ifstr); modelname:= spaste(srcprefix, '.model'); restoredname:= spaste(srcprefix, '.restored'); residualname:= spaste(srcprefix, '.residual'); psffilename:=spaste(srcprefix, '.psf'); printf('IF%d Final imaging for %s\n',spid,fldname[ifld]); imagr.make(modelname); # imagr.weight('natural'); # imagr.makeimage(type='psf',image=psffilename) # bmaj:=0.; bmin:=0.; bpa:=0.; # imagr.fitpsf(psf=psffilename, bmaj=bmaj, bmin=bmin, bpa=bpa) imagr.clean(algorithm='clark', niter=clniter, gain=0.1, threshold='0.0Jy', model=modelname, image=restoredname, residual=residualname); }; } imagr.done(); # # Initialize records to store Stokes info # ivalc:= [=]; qvalc:= [=]; uvalc:= [=]; chicalc:= [=]; # # Iterate through all sources # fldtablename:=spaste(msfilename, '/FIELD'); fldtab:= table(fldtablename); fldname:= fldtab.getcol('NAME'); for (field in fldname) { for (spid in 1:2) { ifstr := sprintf('.if%d',spid); fldprefix:=spaste(prefix, field, ifstr); restoredname:= spaste(fldprefix, '.restored'); # # Measure max abs in (I,Q,U,V) # img1:= image(restoredname); pmax:= array(0.0, 4); for (plane in 1:4) { pregion:= drm.box([1,1,plane], [imsize,imsize,plane]); img1.statistics(statsout=pstats, region=pregion); if (abs(pstats.min) > abs(pstats.max)) { pmax[plane]:= pstats.min; } else { pmax[plane]:= pstats.max; }; }; img1.done(); # ival:= pmax[1]; qval:= pmax[2]; uval:= pmax[3]; ifstr := sprintf('IF%d',spid); ivalc[ifstr][field]:= ival; qvalc[ifstr][field]:= qval; uvalc[ifstr][field]:= uval; } }; # # Print RL phase differences # # angle calibrator as reference # for (spid in 1:2) { ifstr := sprintf('IF%d',spid); qcalc:= qvalc[ifstr][angcal]; ucalc:= uvalc[ifstr][angcal]; chicalc[ifstr]:= atan2 (qcalc, ucalc); } for (field in fieldlist) { print; for (spid in 1:2) { ifstr := sprintf('IF%d',spid); ivala:= ivalc[ifstr][field]; qvala:= qvalc[ifstr][field]; uvala:= uvalc[ifstr][field]; pvala:= sqrt (qvala^2 + uvala^2); chivala:= atan2 (qvala, uvala) - chicalc[ifstr] + absang; chivala:=regang(chivala); print field, " ", ifstr, " I= ", sprintf("%7.3f", ivala), " P= ", sprintf("%7.3f", pvala), " chi= ", sprintf("%7.2f", chivala); } }; # # Successful return value should be T # return T; } killant := function(msf,ants=[]) { # # Use flagger to flag a set of antennas # if( len(ants) < 1 ) { print 'Blank antenna list - no flagging' return } # Assume ms has been created - check with default os tool if( dos.fileexists(file=msf) ) { print 'Found MS file ',msf; } else { fail "MS file does not exist"; } msfilename := msf; kfg := flagger(msfile=msf); kfg.setflagmode('flag'); kfg.setantennas(ants); xf := kfg.query(query='FIELD_ID >= 0'); if ( is_fail(xf) ) { print 'Failure on flagging'; fail; } kfg.done(); print 'Flagged MS file ',msf; # # Successful return value should be T # return T; } quackall := function(msf) { # # Quack first 10s from each scan # # Assume ms has been created - check with default os tool if( dos.fileexists(file=msf) ) { print 'Found MS file ',msf; } else { fail "MS file does not exist"; } msfilename := msf; kfg := flagger(msfile=msf); kfg.setflagmode('flag'); kfg.quack(); kfg.done(); print 'Quacked MS file ',msf; } # # Utility function: atan2(x,y) # Input cos(phi),sin(phi) Output atan in range -180 to +180 # atan2 := function(x, y) { fact:= 180.0 / pi; if ((x > 0) & (y > 0)) angle:= atan(y/x) * fact; if ((x > 0) & (y < 0)) angle:= atan(y/x) * fact + 360; if ((x < 0) & (y > 0)) angle:= 180 - atan(y/abs(x)) * fact; if ((x < 0) & (y < 0)) angle:= 180 + atan(abs(y)/abs(x)) * fact; return angle; }; # # Utility function: regang(x) # regularize angle in degrees to lie from -180 to +180 # regang := function(x) { if ( x > 180.0 ) { while ( x > 180.0 ) x -:= 360.0; } else { if ( x < -180.0 ) { while ( x < -180.0 ) x +:= 360.0; } } return x; }; # # Function to apply gain and opacity correction to ms # vlagaincor := function( par ) { # # Input record parameters: par.msf MS file name (string) # par.tau0 Zenith optical depth (float) # par.band VLA band name (string) # # This function modifies the DATA column of the specified MS. # IT IS RECOMMENDED THAT THE USER MAKE A BACKUP COPY OF THE ORIGINAL MS! # msf := par.msf; tau0 := par.tau0; band := par.band; # # Fill gain coefficients # c := vlagainfill( band ); # # Set the antenna number lookup # anam := spaste(msf, '/ANTENNA'); atab := table(anam); nant := atab.nrows(); anam := atab.getcol("NAME"); printf(' Found %d antennas in table\n',nant); atab.close(); # # Set the frame to be the appropriate telescope # onam := spaste(msf, '/OBSERVATION'); otab := table(onam); tel := otab.getcell("TELESCOPE_NAME",1); printf(' Data is from telescope %s\n',tel); telpos:=dm.observatory(tel); dm.measure(telpos,'WGS84'); # convert to geodetic (why????) dm.doframe(telpos); otab.close(); # # Get source directions # fnam := spaste(msf, '/FIELD'); ftab := table(fnam); nfld := ftab.nrows(); printf(' Found %d fields\n',nfld); fref := ftab.getcolkeyword("REFERENCE_DIR","MEASINFO"); printf(' Reference frame %s\n',fref.Ref); fdirs := ftab.getcol("REFERENCE_DIR"); fields := ftab.getcol("NAME"); ftab.close(); # # Process directions - store in record of directions # dirq := [=]; for ( i in 1:nfld ) { raq := dq.quantity(fdirs[1,,i],'rad'); deq := dq.quantity(fdirs[2,,i],'rad'); dirq[i] := dm.direction(fref.Ref,raq,deq); } # # Make ms table # mst := table(msf, readonly=F); # # Make ms table iterator and iterate # ti := tableiterator(mst, "TIME"); while ( ti.next() ) { st := ti.table(); t := st.getcell("TIME",1); tq := dq.quantity(t,'s'); tymd := dq.time(tq,form='ymd'); print ' Processing timestamp ',tymd; tim := dm.epoch('tai',tymd) ; dm.doframe(tim); # proxy source as first row source for now, get EL for this set ifld := st.getcell("FIELD_ID",1)+1; # (+1 corrects c++ indexing) azel := dm.measure(dirq[ifld],'AZEL'); el := azel.m1; # probably should use getvalue # # Get visibilties - this will be vector [shape=[4 1 nvis] ] # nvis := st.nrows(); vis := st.getcol("DATA"); # # Get AN number vectors (should do something about invalid antennas) # an1v := st.getcol("ANTENNA1")+1; an2v := st.getcol("ANTENNA2")+1; an1 := as_integer( anam[ an1v ]); an2 := as_integer( anam[ an2v ]); # # Construct el numeric vector # elv := array( el.value, nvis ); # # construct opacity corrections = ( 1 - tau/sin(EL) )^(-1) # if ( tau0 > 0.0 ) { oc := opac( tau0, elv ); ocavg := sum( oc )/nvis; } else { oc := array( 1.0, nvis ); ocavg := 1.0; } # # construct gain corrections # gc := gaincv( an1, an2, elv, c ); gcavg := nvis/sum( gc ); # combine these vc := oc/gc; vcavg := sum( vc )/nvis; printf(' Mean correction opac %g , gain %g\n', ocavg, gcavg ); # # make compatible with vis vector # vcv := array(1,4,nvis); for ( i in 1:nvis ) { vcv[1:4,i] := array( vc[i], 4 ); } # # apply to data vector # vis *:= vcv; # # put back in ms # st.putcol("DATA",vis); # NOTE - following 2 lines added 2001-12-18, see 2001-12-18-myers.txt st.flush(); st.done(); printf(' Corrected %s with %d rows at EL = %.2f\n', fields[ifld], nvis, dq.convert(el,'deg').value ); } # # Done # ti.done(); mst.close(); print ' Completed gain and opacity correction\n'; } # # Define the opacity correction # opac := function( tau0, el ) { # # Do atmospheric absorption correction on vector of el # tau(el) = tau0/sin(el) -> oc = 1/[ 1 - tau(el) ] # maximum tau = 1e3 # NOTE - following corrected 2001-12-18, see 2001-12-18-myers.txt tau := tau0/sin( el ); # Cap on high optical depths tau[ tau > 0.999 ] := 0.999; oc := 1.0/( 1.0 - tau ); return oc; } # # Define a gain curve (would index these by antenna eventually) # gaincv := function( an1, an2, el, c) { # # Input vectors (same length) an1,an2,el # Assumes gain curve information is stored in array c[3,nant] # Elevation el in radians # Antenna number an1,an2 elements must correspond to columns in c # Example: Kband AN2 IF1 L (Oct00) # c[1:3,2] := [1.0041e+00, -6.0815e-04, 2.2647e-05 ]; # Note - these coefficients assume ZA in degrees za := 90.0 - el*57.2957795; # Need to pull out appropriate coefficients array pick g1 := c[1,an1] + c[2,an1]*za + c[3,an1]*za*za; g2 := c[1,an2] + c[2,an2]*za + c[3,an2]*za*za; g := sqrt( g1*g2 ); # Need to sort out abnormally low or high gains g[ g < 1.e-4 ] := 1.e-4; g[ g > 1.e4 ] := 1.e4; # done return g; } # # Define function to fill array of gain coefficients # vlagainfill := function( band ) { # fill array as vector if ( band == 'K' | band == 'k' ) { c := [ 1.0645E+00, -2.1800E-03, 1.9091E-06, # AN 1 9.9400E-01, 7.1267E-04, -2.1991E-05, # AN 2 9.0775E-01, 3.4519E-03, -3.2439E-05, # AN 3 8.7106E-01, 4.4819E-03, -3.9123E-05, # AN 4 9.5978E-01, 2.0537E-03, -2.6473E-05, # AN 5 9.6549E-01, 1.7126E-03, -2.1309E-05, # AN 6 9.0808E-01, 3.2659E-03, -2.9129E-05, # AN 7 9.2997E-01, 2.7050E-03, -2.6238E-05, # AN 8 9.4561E-01, 2.4109E-03, -2.7114E-05, # AN 9 9.3171E-01, 2.4068E-03, -2.1307E-05, # AN 10 9.9932E-01, 3.5149E-04, -2.0688E-05, # AN 11 8.4878E-01, 4.8514E-03, -3.9183E-05, # AN 12 9.1357E-01, 2.8711E-03, -2.4172E-05, # AN 13 8.9678E-01, 4.0089E-03, -3.9249E-05, # AN 14 9.8260E-01, 1.4266E-03, -3.0393E-05, # AN 15 8.7698E-01, 4.1427E-03, -3.5156E-05, # AN 16 8.5035E-01, 6.7215E-03, -7.5902E-05, # AN 17 9.9495E-01, 6.8603E-04, -2.4564E-05, # AN 18 9.7246E-01, 1.7074E-03, -2.6990E-05, # AN 19 9.3173E-01, 3.3235E-03, -4.0633E-05, # AN 20 1.0011E+00, 2.3808E-04, -2.2524E-05, # AN 21 9.1979E-01, 3.0976E-03, -3.0121E-05, # AN 22 9.8703E-01, 9.1300E-04, -1.6280E-05, # AN 23 9.9937E-01, 2.3530E-04, -1.4748E-05, # AN 24 8.5658E-01, 4.8396E-03, -4.1113E-05, # AN 25 9.7852E-01, 8.5349E-04, -8.7236E-06, # AN 26 9.2672E-01, 3.0128E-03, -3.1142E-05, # AN 27 9.4566E-01, 2.9426E-03, -4.0201E-05, # AN 28 1.0000E+00, 1.3900E-03, 0.0000E+00, # AN 29 9.4020E-01, 2.4481E-03, -2.9080E-05 ] # AN 30 = avg print ' Filled gain coefficients for VLA at Kband'; } else { if ( band == 'Q' | band == 'q' ) { c := [ 8.5390E-01, 6.9920E-03, -9.7635E-05, # AN 1 na = avg 8.5390E-01, 6.9920E-03, -9.7635E-05, # AN 2 na = avg 8.2699E-01, 8.6734E-03, -1.1053E-04, # AN 3 8.5282E-01, 6.9641E-03, -8.3056E-05, # AN 4 7.8219E-01, 6.6175E-03, -4.7822E-05, # AN 5 8.6940E-01, 7.1348E-03, -1.0007E-04, # AN 6 1.0000E+00, 0.0000E+00, 0.0000E+00, # AN 7 no rx 7.9055E-01, 1.0010E-02, -1.2145E-04, # AN 8 1.0000E+00, 0.0000E+00, 0.0000E+00, # AN 9 no rx 7.9209E-01, 7.7903E-03, -7.4262E-05, # AN 10 8.6475E-01, 3.7786E-03, -6.3711E-05, # AN 11 7.6261E-01, 1.1277E-02, -1.3603E-04, # AN 12 9.5367E-01, 3.4600E-03, -6.8573E-05, # AN 13 8.2851E-01, 9.2871E-03, -1.3003E-04, # AN 14 1.0000E+00, 0.0000E+00, 0.0000E+00, # AN 15 no rx 7.8560E-01, 1.0037E-02, -1.1899E-04, # AN 16 8.5390E-01, 6.9920E-03, -9.7635E-05, # AN 17 na = avg 8.7915E-01, 8.0377E-03, -1.7129E-04, # AN 18 8.7289E-01, 7.1793E-03, -1.0484E-04, # AN 19 7.8001E-01, 9.2160E-03, -9.7229E-05, # AN 20 1.2366E+00, -6.5267E-03, -2.5926E-05, # AN 21 7.5337E-01, 1.1281E-02, -1.3046E-04, # AN 22 8.5390E-01, 6.9920E-03, -9.7635E-05, # AN 23 na = avg 1.0165E+00, -7.3935E-04, -3.9395E-05, # AN 24 6.7155E-01, 1.3864E-02, -1.4794E-04, # AN 25 1.0983E+00, -3.1293E-03, -3.4942E-05, # AN 26 7.9194E-01, 9.6426E-03, -1.1291E-04, # AN 27 9.1300E-01, 5.5659E-03, -9.4216E-05, # AN 28 5.8500E-01, 1.1650E-02, 6.4210E-05, # AN 29 8.5390E-01, 6.9920E-03, -9.7635E-05 ] # AN 30 = avg print ' Filled gain coefficients for VLA at Qband'; } else { if ( band == 'U' | band == 'u' ) { c := [ 1.0213E+00, -1.6223E-04, -1.5142E-05, # AN 1 9.9380E-01, 5.2939E-04, -1.1844E-05, # AN 2 9.5077E-01, 1.7396E-03, -1.5404E-05, # AN 3 9.3945E-01, 1.6803E-03, -9.2158E-06, # AN 4 9.6259E-01, 1.5072E-03, -1.5226E-05, # AN 5 9.7797E-01, 9.2114E-04, -9.6542E-06, # AN 6 9.5072E-01, 1.6461E-03, -1.3873E-05, # AN 7 9.6315E-01, 1.3102E-03, -1.1699E-05, # AN 8 9.8908E-01, 7.4323E-04, -1.2785E-05, # AN 9 9.4594E-01, 1.9071E-03, -1.6902E-05, # AN 10 9.7378E-01, 1.2212E-03, -1.4299E-05, # AN 11 9.2886E-01, 2.0099E-03, -1.4326E-05, # AN 12 9.7188E-01, 1.2060E-03, -1.2958E-05, # AN 13 9.5259E-01, 1.5981E-03, -1.3701E-05, # AN 14 9.9007E-01, 7.6713E-04, -1.5014E-05, # AN 15 9.3265E-01, 2.0155E-03, -1.5228E-05, # AN 16 9.4203E-01, 2.3229E-03, -2.3322E-05, # AN 17 9.8450E-01, 1.3243E-03, -2.9131E-05, # AN 18 9.4425E-01, 2.1415E-03, -2.0617E-05, # AN 19 9.3657E-01, 2.9488E-03, -3.4398E-05, # AN 20 9.9317E-01, 5.8221E-04, -1.4247E-05, # AN 21 9.4996E-01, 1.7295E-03, -1.5152E-05, # AN 22 9.7844E-01, 9.0703E-04, -9.8916E-06, # AN 23 9.7704E-01, 8.9316E-04, -9.0036E-06, # AN 24 9.1715E-01, 2.1045E-03, -1.3433E-05, # AN 25 9.9063E-01, 4.9646E-04, -6.6641E-06, # AN 26 9.5552E-01, 1.6396E-03, -1.5122E-05, # AN 27 9.8287E-01, 8.9881E-04, -1.1829E-05, # AN 28 1.0000E+00, 6.9500E-04, 0.0000E+00, # AN 29 9.6350E-01, 1.3945E-03, -1.5014E-05 ] # AN 30 = avg print ' Filled gain coefficients for VLA at Uband'; } else { c := array( [1,0,0],3,30 ); print ' Filled gain coefficients for VLA to unity'; } } } # reshape as array c::shape := [3,30]; return c; }