#JEH> John Hibbard #JEH> merger.cv.nrao.edu (linux box) #JEH> /aips++/release/aipsinit.csh #JEH> Aug 27 - Nov 11 2001 #JEH> ---------------------- #JEH> First test log from aips++ data reduction #JEH> of NGC 4038 D-array HI spectral line data. #JEH> Week of Aug 27-31 2001 #JEH> Starting testing of spectral line data reduction in aips++ #JEH> NORMAL: Server started: misc (AIPS++ version: 1.6 (build #195)) #JEH> Glish version 2.7. #JEH> #JEH> First, test script 'ngc5921vla.g'. Got data from #JEH> http://www.aoc.nrao.edu/~smyers/aips++/data.html #JEH> Files LINE.UVFITS, CH0.UVFITS - include 'ngc5921vla.g' #JEH> Oh, this seems to execute the script. Save comments to specline1.log #JEH> Map has stripes in ch1-6, 60-63, also lower level stripes in most #JEH> channels. #JEH> Week of September 26-?? #JEH> ----------------------- #JEH> Trying to do with DnC-array data of NGC 4038/9 (AH582, 96JUN03) #JEH> Mode 1A, BWcode=4, 3.125 MHz BW, on-line Hanning #JEH> 128 ch, dv=24.414 kHz ~5.2 km/s #JEH> Flux calibrator =1331+305 (3c286) 15 Jy, 'P', no uvlim #JEH> Phase calibrator=1130-148, dr=8.44deg, 6.2 Jy, 'P', no uvlim #JEH> 30sec integrations on calibrators #JEH> 60sec integrations on source #JEH> #JEH> Notes from AIPS data reduction: #JEH> refant=24 #JEH> bl 4-10 bad #JEH> some 1st scans #JEH> ant 10 until 0/1 12 0 #JEH> ant 17 has 18deg jump #JEH> ant 18 flagged by on-line system (W1=shadowed) #JEH> 2 kJy source on line data, bl 8-23 #JEH> self-cal made a big difference (3 phase only, 1 ph+gain) #JEH> ch 0 noise went from 1 mJy/b to 0.4 mJy/b (2300 CC) #JEH> unlin ch 5-23+105-115 #JEH> R+1: 47"x31", cleaned to 2.5 mJy/b. 1sig=1.1 mJy/b #JEH> #JEH> Had analysts run fillm. File=AH582.960903.LINE (118 MB) #JEH> #JEH> ======================================================== #JEH> USING n5921.g script and Ch5,6 of Getting Results v.Sept7,2001 #JEH> => nope, ch5,6 are useless. What am I supposed to use? #JEH> Can't find any scripts that include calibration. #JEH> Getting Started and Getting Results are useless. #JEH> Using John Bensons guide instead. #JEH> 9/27: appendix C of Getting Results (9/27/01 version) #JEH> "VLA reduction in AIPS++" #JEH> Also using some of Greg and Steves test logs as a guide. #JEH> Start of glish: - include 'ms.g' - fitsfilename :='/export/data_1/jhibbard/DATA/n4038/AH582/AH582.960903.LINE'; - myms:=fitstoms(fitsfile=fitsfilename, msfile='n4038D.ms'); - myms.summary(verbose=T); #GUI> ================================================================ #GUI> Equivalent Gui Commands: #GUI> in Tool Manager window, select general.ms.ms #GUI> select 'fitstoms' constructor #GUI> toolname=myms #GUI> msfile=n4038D.ms #GUI> filenameFrom catalog #GUI> (make sure to specify Type as "FITS") #GUI> #GUI> readonly TRUE #GUI> lock FALSE #GUI> #GUI> (didnt take too long) #GUI> (Tool Manager now shows the functions for myms) #GUI> select 'summary' function #GUI> verbose TRUE #GUI> #GUI> ================================================================ #JEH>Logger messages after creating myms: #LOG> ================================================================ #LOG>Starting filecatalog GUI. Select an entry #LOG>Execution will construct new tool myms #LOG>Converting FITS file '/export/data_1/jhibbard/DATA/n4038/AH582/AH582.960603.LINE' to MeasurementSet 'n4038D.ms' #LOG>Reading and writing 77213 visibility groups #LOG>Found binary table of type AIPS AN following data #LOG>Found binary table of type AIPS NX following data #LOG>Skipping table, duplicate or unrecognized type: AIPS NX #LOG>Found binary table of type AIPS SU following data #LOG>Found binary table of type AIPS FQ following data #LOG>Found binary table of type AIPS CL following data #LOG>Skipping table, duplicate or unrecognized type: AIPS CL #LOG>Found binary table of type AIPS TY following data #LOG>Skipping table, duplicate or unrecognized type: AIPS TY #LOG> #LOG>Logger messages after summary function: #LOG> #LOG>Flushing MS to disk #LOG> #LOG> MeasurementSet Name: n4038D.ms MS Version 2 #LOG> #LOG> Observer: AH582 Project: #LOG>Observation: VLA #LOG>Fields: 3 #LOG> ID Name Right Ascension Declination Type #LOG> 1 1331+305 13:31:08.29 +30.30.32.96 J2000 #LOG> 2 1130-148 11:30:07.05 -14.49.27.39 J2000 #LOG> 3 NGC4038 12:01:52.00 -18.55.40.00 J2000 #LOG>Data records: 77213 Total integration time = 12450 seconds #LOG> Observed from 03-Jun-1996/00:56:00 to 03-Jun-1996/04:23:30 #LOG>Data descriptions: 1 (1 spectral windows and 1 polarization setups) #LOG> ID Ref.Freq #Chans Resolution TotalBW Correlations #LOG> 1 1412.63MHz 127 24.4141kHz 3125 kHz RR #LOG>Feeds: 29: printing first row only #LOG> Antenna Spectral Window # Receptors Polarizations #LOG> 1 -1 2 [ R, L] #LOG>Antennas: 29: #LOG> ID Name Station Diam. Long. Lat. #LOG> 1 1 VLA:N18 25.0 m -107.37.12.0 +33.54.58.3 #LOG> 2 2 VLA:W7 25.0 m -107.37.18.4 +33.53.54.8 #LOG> 3 3 VLA:W2 25.0 m -107.37.07.4 +33.54.00.9 #LOG> 4 4 VLA:E8 25.0 m -107.36.48.9 +33.53.55.1 #LOG> 5 5 VLA:W9 25.0 m -107.37.25.1 +33.53.51.0 #LOG> 6 6 VLA:E1 25.0 m -107.37.05.7 +33.53.59.2 #LOG> 7 7 VLA:E2 25.0 m -107.37.04.4 +33.54.01.1 #LOG> 8 8 VLA:W8 25.0 m -107.37.21.6 +33.53.53.0 #LOG> 9 9 VLA:N14 25.0 m -107.37.09.9 +33.54.38.5 #LOG> 10 10 VLA:W3 25.0 m -107.37.08.9 +33.54.00.1 #LOG> 11 11 VLA:N4 25.0 m -107.37.06.5 +33.54.06.1 #LOG> 12 12 VLA:W6 25.0 m -107.37.15.6 +33.53.56.4 #LOG> 13 13 VLA:W4 25.0 m -107.37.10.8 +33.53.59.1 #LOG> 14 14 VLA:N6 25.0 m -107.37.06.9 +33.54.10.3 #LOG> 15 15 VLA:N16 25.0 m -107.37.10.9 +33.54.48.0 #LOG> 16 16 VLA:E4 25.0 m -107.37.00.8 +33.53.59.7 #LOG> 17 17 VLA:E3 25.0 m -107.37.02.8 +33.54.00.5 #LOG> 18 18 VLA:W1 25.0 m -107.37.05.9 +33.54.00.5 #LOG> 19 19 VLA:E6 25.0 m -107.36.55.6 +33.53.57.7 #LOG> 20 20 VLA:N10 25.0 m -107.37.08.2 +33.54.22.4 #LOG> 21 21 VLA:E9 25.0 m -107.36.45.1 +33.53.53.6 #LOG> 22 22 VLA:E7 25.0 m -107.36.52.4 +33.53.56.5 #LOG> 23 23 VLA:N12 25.0 m -107.37.09.0 +33.54.30.0 #LOG> 24 24 VLA:N2 25.0 m -107.37.06.2 +33.54.03.5 #LOG> 25 25 VLA:OUT 25.0 m -107.37.06.0 +33.54.01.7 #LOG> 26 26 VLA:E5 25.0 m -107.36.58.4 +33.53.58.8 #LOG> 27 27 VLA:N8 25.0 m -107.37.07.5 +33.54.15.8 #LOG> 28 28 VLA:W5 25.0 m -107.37.13.0 +33.53.57.8 #LOG> 29 29 VLA:VPT 25.0 m -107.37.06.0 +33.54.01.7 #LOG> #LOG> #LOG>Tables: #LOG> MAIN 77213 rows #LOG> ANTENNA 29 rows #LOG> DATA_DESCRIPTION 1 row #LOG> DOPPLER #LOG> FEED 29 rows #LOG> FIELD 3 rows #LOG> FLAG_CMD #LOG> FREQ_OFFSET #LOG> HISTORY 29 rows #LOG> OBSERVATION 1 row #LOG> POINTING 319 rows #LOG> POLARIZATION 1 row #LOG> PROCESSOR #LOG> SOURCE (see FIELD) #LOG> SPECTRAL_WINDOW 1 row #LOG> STATE #LOG> SYSCAL #LOG> WEATHER #LOG> #LOG>=================================================================== #JEH> Wanted to get a listr-like output, so tried: - myms.lister(startime="03-Jun-1996/00:56:00" , stoptime="03-Jun-1996/04:23:30" ); #JEH> (times pasted in from summary output) #GUI> ================================================================ #GUI> Equivalent Gui Commands: #GUI> select 'lister' function #GUI> startime 03-Jun-1996/00:56:00 #GUI> stoptime 03-Jun-1996/04:23:30 #GUI> #GUI> ================================================================ #JEH> well, that didn't seem to work. It produced the exact same #JEH> output as myms.summary, then went away for a long time then #JEH> aborted with: #BUG> - LocalExec::SetStatus: abnormal child termination for /aips++/stable/linux/bin/ms #BUG> Server 'ms' has failed unexpectedly! #BUG> You will need to create the relevant tool again. #BUG> If that causes unexpected behavior, please restart AIPS++ #BUG> Please submit a bug-report using bug() if you can reproduce the problem. #JEH> I looked at Tools in use, and myms is still there, so I don't #JEH> know what the complaint it. #JEH> Selected myms from Tools in use, , then repeated summary to #JEH> make sure it is there. Went away for a long time. Too long. >5min. #JEH> So, running abort, delete myms, and will recreate. #JEH> Hmmm - doesn't kill: #LOG> Unable to find server for tool myms since it does not have an id() function #JEH> Went to "Tools in use", selected myms and #JEH> Recreate myms. #JEH> The next part of the script tells me to flag using #JEH> flagger. However, you would have to know before hand #JEH> what needs flagging, and I don't, so lets use msplot instead #JEH> oops - this is all the data! I need to do the equivalent of #JEH> AVSPC. I think that would be myms.selectchannel, but there is #JEH> no place for an output. Perhaps I need to create a NEW ms tool, #JEH> and use that. I'll try. #JEH> Hmm....nothing happening. Not getting the filler progress #JEH> meter. Something is screwy. Can't get anything to come back. #JEH> Go back to Tool manager, and can create an ms tool, but nothing #JEH> responds. Guess I have to abort. # #JEH> Ah, I see that Steve Myers already bugged this in his 2001.05.02 #JEH> log file. Related to the earlier ms.lister bug. Ed says its b/c #JEH> lister can only take list time periods. #JEH> So I'll leave it alone, abort, and restart. include 'ms.g'; myms:=ms(filename="/export/data_1/jhibbard/aips++/n4038D.ms" , readonly=T, lock=F); ok:=myms.selectchannel(nchan=1, start=15, width=96, inc=1); #JEH> Comes back immediately. So really doesn't do anything, does it? myms.done(); #JEH> George says we need to flag the autocorrelation data. - include 'flagger.g' - fg:=flagger(msfile="/export/data_1/jhibbard/aips++/n4038D.ms" ); - fg.auto(); - fg.done(); #GUI> ================================================================ #GUI> synthesis.flagger.flagger #GUI> select constructor flagger #GUI> toolname myflagger #GUI> msfile=/export/data_1/jhibbard/aips++/n4038D.ms #GUI> #GUI> select function auto #GUI> trial false #GUI> #GUI> |done| #GUI> ================================================================ #JEH> Had some problems with Tk/Tcl and RH7.1; P.Murphy ran a script #JEH> he got from Athol, so I can now run weekly. So I am not in #JEH> Weekly (build#279). That's not a problem is it? #JEH> Lets take a look at the data. - include 'msplot.g' - mymsplot:=msplot('n4038D.ms', edit=T, flagfile="n4038Dflag01.flag"); #GUI> ================================================================ #GUI> general.ms.msplot #GUI> toolname mymsplot #GUI> msfile n4038D.ms #GUI> edit TRUE #GUI> flagfile n4038Dflag01.flag #GUI> #GUI> ================================================================ #LOG> Starting server ms #LOG> Server started: /aips++/weekly/linux/bin/ms (AIPS++ version: 1.6 (build #279)) #LOG> Opened the measurement set called n4038D.ms for read/write access. #LOG> #LOG> MeasurementSet Name: n4038D.ms MS Version 2 #LOG> #LOG> Observer: AH582 Project: #LOG> Observation: VLA(29 antennas) #LOG> Fields: 3 #LOG> Name: 1331+305 RA: 13:31:08.29 Dec: +30.30.32.96 (J2000 ) #LOG> Name: 1130-148 RA: 11:30:07.05 Dec: -14.49.27.39 (J2000 ) #LOG> Name: NGC4038 RA: 12:01:52.00 Dec: -18.55.40.00 (J2000 ) #LOG> Data records: 77213 Total integration time = 12450 seconds #LOG> Observed from 03-Jun-1996/00:56:00 to 03-Jun-1996/04:23:30 #LOG> Data descriptions: 1 (1 spectral windows and 1 polarization setups) #LOG> ID Ref.Freq #Chans Resolution TotalBW Correlations #LOG> 1 1412.63MHz 127 24.4141kHz 3125 kHz RR #LOG> Antennas: 29 #LOG> ID= 1-6: 1=VLA:N18, 2=VLA:W7, 3=VLA:W2, 4=VLA:E8, 5=VLA:W9, 6=VLA:E1, #LOG> ID= 7-12: 7=VLA:E2, 8=VLA:W8, 9=VLA:N14, 10=VLA:W3, 11=VLA:N4, 12=VLA:W6, #LOG> ID= 13-17: 13=VLA:W4, 14=VLA:N6, 15=VLA:N16, 16=VLA:E4, 17=VLA:E3, #LOG> ID= 18-22: 18=VLA:W1, 19=VLA:E6, 20=VLA:N10, 21=VLA:E9, 22=VLA:E7, #LOG> ID= 23-27: 23=VLA:N12, 24=VLA:N2, 25=VLA:OUT, 26=VLA:E5, 27=VLA:N8, #LOG> ID= 28-29: 28=VLA:W5, 29=VLA:VPT #LOG> #LOG> #LOG> Tables(rows): (-1 = table absent) #LOG> MAIN(77213) #LOG> ANTENNA(29) DATA_DESCRIPTION(1) DOPPLER(-1) FEED(29) FIELD(3) #LOG> FLAG_CMD(0) FREQ_OFFSET(-1) HISTORY(29) OBSERVATION(1) POINTING(319) #LOG> POLARIZATION(1) PROCESSOR(0) SOURCE(-1) SPECTRAL_WINDOW(1) STATE(0) #LOG> SYSCAL(-1) WEATHER(-1) #LOG> #LOG> Copying the flag columns to a table called n4038Dflag01.flag #LOG> This may take some time... #LOG> Reading MS to get ranges for selected data: this may take a while... #LOG> ...found all ranges #JEH> Works entirely from GUI: #GUI> Show: Plot X vs Y, #GUI> Select: X axis=uvdist, Yaxis=data #GUI> Data selection: Field Id=1 #GUI> #LOG> ********** Plot Observed data amplitude versus UV Distance ********** #LOG> Selecting a subset of the measurement set. #LOG> This may take a while, particularly for a complicated selection. #LOG> ...subset has 14040 rows. #LOG> Channel selection: #chan=127, start=1, width=1, incr=1 #LOG> Polarization selection: [RR] #LOG> Selection initialized ok #LOG> Polarization selection: [RR] #LOG> Channel selection: #chan=127, start=1, width=1, incr=1 #LOG> Reading by time #LOG> Plotted 1783080 points on 1 pages #JEH> I guess the ms.setchannel didn't keep. #JEH> Most data around zero; one point at > 100 Jy #GUI> do box around high point and press #JEH> Doesn't rescale and it really should. #JEH> Ooops - didn't press the Flag all Channels correlations. #GUI> select Flag all channels, correlations. #JEH> OK, now it rescaled. #JEH> Repeat for phase cal: #GUI> Data selection: Field Id=2 #GUI> #JEH> Ok, Repeat for n4038: #GUI> Data selection: Field Id=3 #GUI> #LOG> Plotted 6527038 points on 1 pages #JEH> Most data around zero; one point at > 100 Jy #GUI> do box around high point and press #JEH> Hmmm. This time it didn't rescale. Wonder why. #GUI> I was trying to zoom in on the plot by pressing various mouse #GUI> buttons, and got in a weird state that I can't seem to get out #GUI> of. Time to read the documentation. #JEH> Really this is ridiculous. Only need to see an average of #JEH> the channels here. #JEH> Wait, under the Spectral selection: #GUI> List of spectral windows [1] #GUI> Number of channels 1 #GUI> Start channel 11 #GUI> Width of channel 96 #GUI> Channel increment 1 #GUI> ***DONT FORGET TO PRESS RETURN AFTER ENTERING A VALUE!*** #LOG> Selection initialized ok #LOG> Polarization selection: [RR] #LOG> Channel selection: #chan=1, start=11, width=96, incr=1 #LOG> Reading by time #LOG> Plotted 51389 points on 1 pages #JEH> MUCH MUCH better. #JEH> Time to go; more tomorrow. #JEH> 9/28/01 #JEH> Not sure what flags got applied when I was doing the zoom/unzoom #JEH> stuff, so pressed "Revert". Never came back: #LOG> Caught an exception! Event type=run exception=Array::operator()(b,e,i) - b,e or i incorrectly specified #LOG> Caught an exception! Event type=run exception=Array::operator()(b,e,i) - b,e or i incorrectly specified #LOG> Caught an exception! Event type=run exception=Array::operator()(b,e,i) - b,e or i incorrectly specified #LOG> Caught an exception! Event type=run exception=Array::operator()(b,e,i) - b,e or i incorrectly specified #LOG> Caught an exception! Event type=run exception=Array::operator()(b,e,i) - b,e or i incorrectly specified #BUG> Had to cntr-c out of aips++ #BUG> Submitted Bug AOCso03034 on above #BUG> Also submitted enhancement to suggest keeping flagging commands #BUG> instead of entering 0's & 1's, and to ask what a row is. #BUG> enchancment AOCso03035 #JEH> Lets start again - include 'msplot.g' - mymsplot:=msplot('n4038D.ms', edit=T, flagfile="n4038Dflag01.flag"); #GUI> Show: Plot X vs Y, #GUI> Select: X axis=uvdist, Yaxis=data #JEH> Spectral selection: #GUI> List of spectral windows [1] #GUI> Number of channels 1 #GUI> Start channel 11 #GUI> Width of channel 96 #GUI> Channel increment 1 #GUI> #LOG> Polarization selection: [RR] #LOG> Channel selection: #chan=1, start=11, width=96, incr=1 #LOG> Reading by time #LOG> Plotted 77213 points on 1 pages #JEH> Thats strange; no outliers now. Maybe averaging hides them? #JEH> Or did the flags keep? I deleted the earlier flag file #JEH> Spectral selection: #GUI> List of spectral windows [1] #GUI> Number of channels 127 #GUI> Start channel 1 #GUI> Width of channel 1 #GUI> Channel increment 1 #GUI> #JEH> Ah, there are the outliers - so at the edges of the band. #JEH> Marked them and flagged them. Takes forever. #JEH> goddammit. Forgot to set the all channels, all correlations #JEH> flags again. Had to repeat. #JEH> OK, good things to look at in msplot: #JEH> X=uvdist, Y=data #JEH> X=time, Y=data #JEH> Saw some low data at various times, and wanted to see who they #JEH> were, so typed locate, but nothing came back: #LOG> Region is [4180.50146 0.263511658] to [7477.24561 0.177205995] #LOG> Channel selection: #chan=127, start=1, width=1, incr=1 #LOG> Polarization selection: [RR] #LOG> Selection initialized ok #LOG> Polarization selection: [RR] #LOG> Channel selection: #chan=1, start=11, width=96, incr=1 #LOG> Reading by time #BUG> Submitted bug AOCso03036 #JEH> Also noticed that the "Flag all channels" box becomes inactive. #JEH> I want any flagging to flag all channels, but I don't want to #JEH> have to plot all points. I think this is the scope problem #JEH> others have talked about, but I'll put it anyways. #BUG> submitted Bug AOCso03037 #JEH> X=time, Y=data, dataselect fieldid, iterate on antennae 1 #JEH> Hmm. This will be much more useful after preliminary #JEH> calibration. One last thing; lets try the tvflag type plot: #JEH> Show: display data as an image #JEH> Xaxis=time, Y-axis=ifr_number, Image pixels=observed #JEH> plots interferometer_number in km vs Time in km?????? #JEH> Cant tell anything from this. No idea what is really being #JEH> plotted. Dismiss wouldn't get rid of it, so had to "x" it. #BUG: submitted Bug AOCso03038 #JEH> Ok, enough of this. Lets try priliminary calibration. #JEH> msplot.done(); - include 'imager.g'; - im:=imager('n4038D.ms'); #JEH> Following along in the n5921 script, I see its time to #JEH> set the flux. First there is a call to setjy with defaults. #JEH> Does this set all fluxes to 1 Jy? Bensons notes say we have #JEH> to do this, so lets try it. - im.setjy(); #LOG> Starting imager::setjy #LOG> Selecting data #LOG> Performing selection on MeasurementSet #LOG> Selecting on field and spectral window ids #LOG> No selection string given #LOG> By selection 77213 rows are reduced to 14040 #LOG> Selecting 127 channels, starting at visibility channel 1 stepped by 1 #LOG> Fourier transforming: replacing MODEL_DATA column #LOG> Fourier transforms will use image centers as tangent points #LOG> made generic SkyEquation #LOG> 1331+305 spwid= 1 [I=14.76, Q=0, U=0, V=0] Jy, (Perley-Taylor 95) #LOG> Selecting data #LOG> Performing selection on MeasurementSet #LOG> Selecting on field and spectral window ids #LOG> No selection string given #LOG> By selection 77213 rows are reduced to 11779 #LOG> Selecting 127 channels, starting at visibility channel 1 stepped by 1 #LOG> Fourier transforming: replacing MODEL_DATA column #LOG> made generic SkyEquation #LOG> 1130-148 spwid= 1 [I=1, Q=0, U=0, V=0] Jy, (default) #LOG> Selecting data #LOG> Performing selection on MeasurementSet #LOG> Selecting on field and spectral window ids #LOG> No selection string given #LOG> By selection 77213 rows are reduced to 51394 #LOG> Selecting 127 channels, starting at visibility channel 1 stepped by 1 #LOG> Fourier transforming: replacing MODEL_DATA column #LOG> made generic SkyEquation #LOG> NGC4038 spwid= 1 [I=1, Q=0, U=0, V=0] Jy, (default) #LOG> Finished imager::setjy #LOG> 45.6 real 26.01 user 9.63 system #JEH> So it looks like this call to setjy actually calculates the #JEH> flux of the primary, and sets everything else to 1 Jy. Good. #JEH> However, in aips setjy found S=14.8929 Jy for 3c286, so I'll #JEH> set this instead of the 14.76 set in the above: #JEH> Primary flux call is 1331+305 (3c286), fieldid 1 from summary, - im.setjy(fieldid=1,fluxdensity=14.8929); - im.done(); #LOG> Starting imager::setjy #LOG> Selecting data #LOG> Performing selection on MeasurementSet #LOG> Selecting on field and spectral window ids #LOG> No selection string given #LOG> By selection 77213 rows are reduced to 14040 #LOG> Selecting 127 channels, starting at visibility channel 1 stepped by 1 #LOG> Fourier transforming: replacing MODEL_DATA column #LOG> made generic SkyEquation #LOG> 1331+305 spwid= 1 [I=14.89, Q=0, U=0, V=0] Jy, (user-specified) #LOG> Finished imager::setjy #LOG> 7.64 real 4.69 user 1.77 system #JEH> Continuing on 10/01/01 #JEH> Look to make sure above fluxes got into the proper table: #GUI> #GUI> select n4038D.ms under File Name, click on #GUI> Brings up Table Browser. #GUI> Slide over to the MODEL_DATA colum and click on any row #GUI> This will give you the flux value, but doesn't tell you #GUI> who it is for. We want to see the flux value for #GUI> field_id=[1,2,3], so I'll try to do this using the #GUI> Table->Select pulldown menu #GUI> I do this, selecting MODEL_DATA for FIELD_ID==1, and #GUI> it sends the output to some file that I can't view via dc. #JEH> From viewing the tester logs, it seems as if this should #JEH> be possible via some combination of ms.select and table #JEH> tools, but I can't quite make sense of it. #JEH> Apparently, it is possible to look at this with msplot #JEH> But this is inelegant and non-exact. #JEH> Tried this: - include 'table.g' - t:=table('n4038D.ms'); - fields:=t.getcol("FIELD_ID"); - model:=t.getcol("MODEL_DATA"); - shape(model); [1 127 77213] - s0:=model[1,1,fields==0]; - s1:=model[1,1,fields==1]; - s2:=model[1,1,fields==2]; - s3:=model[1,1,fields==3]; - print s0[1] 14.8929005+0i - print s1[1] 1-5.15758192e-14i - print s2[1] 1+1.01126491e-15i - print s3[1] : index (= 1) out of range, array length = 0 - t.done(); #-------------------------------------- #JEH> George wrote back and said the following: #If you want to avoid dealing with the zero-counting, you can use the #ms tool to access the Measurement Set columns as follows. - include 'ms.g' - myms:=ms('n4038D.ms') - myms.summary(verbose=T) # to get list of fields, antennas, etc. # Get model column for FIELD_ID=1, print out unique entries: - myms.selectinit() # In next line note [] syntax (items is a record, not a string) - myms.select(items=[field_id=1]); - print unique(myms.getdata("MODEL_DATA").model_data); 14.8929005+0i # myms.getdata(..) returns a record with a field called model_data, # and the unique function will limit the values printed to unique ones # Get model column for FIELD_ID=2, print out unique entries: - myms.selectinit() - myms.select(items=[field_id=2]); - print unique(myms.getdata("FIELD_ID MODEL_DATA").field_id); 2 - print unique(myms.getdata("FIELD_ID MODEL_DATA").model_data); # [1-3.91456008e-12i 1-3.91449286e-12i 1-3.91442521e-12i 1-3.91435755e-12i 1-3.9142899e-12i 1-3.91422224e-12i 1-3.91415459e-12i 1-3.91408737e-12i 1-3.91401972e-12i 1-3.91395206e-12i 1-3.91388441e-12i 1-3.91381675e-12i 1-3.9137491e-12i 1-3.91368188e-12i 1-3.91361422e-12i 1-3.91354657e-12i 1-3.91347892e-12i 1-3.91341126e-12i 1-3.91334404e-12i 1-3.91327639e-12i 1-3.91320873e-12i 1-3.91314108e-12i 1-3.91307342e-12i 1-3.91300577e-12i 1-3.91293855e-12i 1-3.9128709e-12i 1-3.91280324e-12i 1-3.91273559e-12i 1-3.91266793e-12i 1-3.91260028e-12i 1-3.91253306e-12i 1-3.9124654e-12i 1-3.91239775e-12i 1-3.9123301e-12i 1-3.91226244e-12i 1-3.91219479e-12i 1-3.91212757e-12i 1-3.91205991e-12i 1-3.91199226e-12i 1-3.9119246e-12i 1-3.91185695e-12i 1-3.9117893e-12i 1-3.91172207e-12i 1-3.91165442e-12i 1-3.91158677e-12i 1-3.91151911e-12i 1-3.91145146e-12i 1-3.9113838e-12i 1-3.91131658e-12i 1-3.91124893e-12i 1-3.91118127e-12i 1-3.91111362e-12i 1-3.91104597e-12i 1-3.91097831e-12i 1-3.91091109e-12i 1-3.91084344e-12i 1-3.91077578e-12i 1-3.91070813e-12i 1-3.91064047e-12i 1-3.91057282e-12i 1-3.9105056e-12i 1-3.91043795e-12i 1-3.91037029e-12i 1-3.91030264e-12i 1-3.91023498e-12i 1-3.91016733e-12i 1-3.91010011e-12i 1-3.91003245e-12i 1-3.9099648e-12i 1-3.90989715e-12i 1-3.90982949e-12i 1-3.90976184e-12i 1-3.90969462e-12i ... # In this case, both the FIELD_ID and MODEL_DATA columns have # been extracted into the record returned by getdata. The record fields # assume the names of the columns extracted. You need to be careful # using getdata because if you ask for too many columns from a large MS, # you will fill up the memory on your machine! #JEH> Here is what I cobbled together, combining both my and Georges #JEH> approaches: - include 'ms.g' - myms:=ms('n4038D.ms') - myms.select(items=[field_id=1]); - s1:=myms.getdata("MODEL_DATA").model_data - shape(s1); [1 127 14040] - print s1[1] 14.8929005+0i - myms.selectinit(); - myms.select(items=[field_id=2]); - s2:=myms.getdata("MODEL_DATA").model_data - shape(s2); [1 127 11779] - print s2[1] 1-5.15758192e-14i - myms.selectinit(); - myms.select(items=[field_id=3]); - s3:=myms.getdata("MODEL_DATA").model_data - shape(s3); [1 127 51394] - print s3[1] 1+1.01126491e-15i - myms.done(); #-------------------------------------- #JEH> So fields here and in summary do not match. Bugged it. #BUG> AIPS++ defect AOCso03046: FIELD_ID for table not same as in ms.summary #BUG> Also enhancement request to list FAQ's in front page: #BUG> AIPS++ defect AOCso03047: should be a summary of all questions in FAQ #JEH> Enough for today. #JEH> 10/02/01. Onward. #JEH> Just following along in Georges' script. #JEH> Says it is calculating G Jones "using psuedo continuum" #JEH> but I have no idea what that means. Lets hope its #JEH> not calculating a separate Jones matrix for each channel. #JEH> Emailed George to request more information. #JEH> From aips, there are no uvlimits for the calibrators, and #JEH> we used refant=24 - include 'calibrater.g'; - cal:=calibrater('n4038D.ms'); #JEH> First, solve for complex gains. The example is written for #JEH> 5min timescale, but I think we want scan averages. Looking #JEH> at my AIPS lister output, I observe 1331+305 for 13min=780sec, #JEH> so I want a long average. Bug this. #BUG> AIPS++ defect AOCso03048: hard wired numbers need to be explained #JEH> NOTE! THE MATRIX BETTER NOT CROSS FIELD_ID BOUNDARIES!! #JEH> Also, no definition of preavg in on-line documentation. Bug it. #BUG> AIPS++ defect AOCso03049: no description of preavg - cal.setdata(msselect='FIELD_ID<=2'); - cal.setsolve(type='G',t=900.0,preavg=900.0,refant=24,table='n4038D.gcal'); - cal.solve(); LOG> Starting calibrater::solve LOG> Initializing solvable electronic gain terms (G-matrix) LOG> For interval of 900 seconds, found 8 slots LOG> Solving for G LOG> G Jones Slot=1, 1331+305, spw=1: 03-Jun-1996/00:45:00 to 03-Jun-1996/00:59:59 LOG> G Jones Initial fit per unit weight = 1.1215 Jy, sum of weights = 11224 LOG> G Jones Final fit per unit weight = 0.0131133 Jy after 8 iterations LOG> G Jones Slot=2, 1331+305, spw=1: 03-Jun-1996/01:00:00 to 03-Jun-1996/01:14:59 LOG> G Jones Initial fit per unit weight = 1.09631 Jy, sum of weights = 25268 LOG> G Jones Final fit per unit weight = 0.0128813 Jy after 8 iterations LOG> G Jones Slot=3, 1331+305, spw=1: 03-Jun-1996/04:15:00 to 03-Jun-1996/04:29:59 LOG> G Jones Initial fit per unit weight = 1.09207 Jy, sum of weights = 19656 LOG> G Jones Final fit per unit weight = 0.0101289 Jy after 8 iterations LOG> G Jones Slot=4, 1130-148, spw=1: 03-Jun-1996/01:00:00 to 03-Jun-1996/01:14:59 LOG> G Jones Initial fit per unit weight = 0.271481 Jy, sum of weights = 5616 LOG> G Jones Final fit per unit weight = 0.0061215 Jy after 4 iterations LOG> G Jones Slot=5, 1130-148, spw=1: 03-Jun-1996/01:45:00 to 03-Jun-1996/01:59:59 LOG> G Jones Initial fit per unit weight = 0.413641 Jy, sum of weights = 11128 LOG> G Jones Final fit per unit weight = 0.00581371 Jy after 8 iterations LOG> G Jones Slot=6, 1130-148, spw=1: 03-Jun-1996/02:30:00 to 03-Jun-1996/02:44:59 LOG> G Jones Initial fit per unit weight = 0.40201 Jy, sum of weights = 10716 LOG> G Jones Final fit per unit weight = 0.00554595 Jy after 8 iterations LOG> G Jones Slot=7, 1130-148, spw=1: 03-Jun-1996/03:15:00 to 03-Jun-1996/03:29:59 LOG> G Jones Initial fit per unit weight = 0.397881 Jy, sum of weights = 9828 LOG> G Jones Final fit per unit weight = 0.00521199 Jy after 8 iterations LOG> G Jones Slot=8, 1130-148, spw=1: 03-Jun-1996/04:00:00 to 03-Jun-1996/04:14:59 LOG> G Jones Initial fit per unit weight = 0.393614 Jy, sum of weights = 9828 LOG> G Jones Final fit per unit weight = 0.00492674 Jy after 8 iterations LOG> Storing G matrix in table n4038D.gcal LOG> Finished calibrater::solve LOG> 30.27 real 21.25 user 5.23 system JEH> Well, that is bullshit. I don't want it to start at 03-Jun-1996/00:45:00; JEH> I want it to start at the start of the observation!! JEH> Maybe this what preavg is? I'll repeat, using preavg of 300s. JEH> Nope; results in same time slots. But initial fits per unit JEH> weight change. Try preavg=0. Nope. JEH> I can't see a way to get what I want. bugit. BUG> AIPS++ defect AOCso03050: time slots are ignorant JEH> I have *NO* idea what this is doing in the spectral dimension. JEH> The defaults in calibrater.setdata are nchan=1,start=1,step=1 JEH> but the log did not plot results for each of the 127 channels. JEH> Call George. No answer. Call Athol. No answer. Email Athol. - cal.plotcal(tablename='n4038D.gcal') JEH> Can see that ant.14 has so very descrepant values. Leave it for now. JEH> First scan on 1130-148 has high gains - this had only 1.5 min on JEH> the calibrator, so is a lower S/N (see weights in above LOG). JEH> Transfer the flux density scale from 1331+305 to 1130-148: - cal.fluxscale(tablein='n4038D.gcal', tableout='n4038D.gcalout', reference='1331+305', transfer=['1130-148']); - cal.plotcal(tablename='n4038D.gcalout') LOG> Set the flux density scale using reference calibrators LOG> Flux density for 1130-148 (spw=1) is: 2.171 +/- 0.002 Jy JEH> AIPS initially found S=5.44449 +/- 0.0705 (5.37044 +/- 0.0274 JEH> after flagging). Quite a difference. Nominal value at L-band JEH> is 6.2 Jy. So what up with this? JEH> Lets proceed through the script to see how bad it screws things up. JEH> Next step is to calculate the bandpass. I used 1331+305 only as BPcal. - cal.setdata(msselect='FIELD_ID==1'); - cal.setapply(type='G',t=0.0,table='n4038D.gcalout'); - cal.setsolve(type='B',t=86400.0,preavg=86400.0,refant=24,table='n4038D.bcal'); JEH> Aborted: LOG> messages in log window: Starting calibrater::solve Initializing solvable bandpass terms (B-matrix) For interval of 86400 seconds, found 1 slots Loading G table from n4038D.gcalout Solving for G G Jones Slot=1, 1331+305, spw=1: 03-Jun-1996/00:45:00 to 03-Jun-1996/00:59:59 G Jones Initial fit per unit weight = 3.52002 Jy, sum of weights = 29472 G Jones Final fit per unit weight = 3.52002 Jy after 1 iterations G Jones Slot=2, 1331+305, spw=1: 03-Jun-1996/01:00:00 to 03-Jun-1996/01:14:59 G Jones Initial fit per unit weight = 0.569396 Jy, sum of weights = 7020 G Jones Final fit per unit weight = 0.0136284 Jy after 4 iterations G Jones Slot=3, 1130-148, spw=1: 03-Jun-1996/01:45:00 to 03-Jun-1996/01:59:59 G Jones Initial fit per unit weight = 0.0114765 Jy, sum of weights = 12636 G Jones Final fit per unit weight = 0.0114765 Jy after 1 iterations G Jones Slot=4, 1130-148, spw=1: 03-Jun-1996/02:30:00 to 03-Jun-1996/02:44:59 G Jones Initial fit per unit weight = 0.0149179 Jy, sum of weights = 7020 G Jones Final fit per unit weight = 0.0149179 Jy after 1 iterations Solving for B B Jones Slot=1, 1331+305, spw=1: 03-Jun-1996/24:00:00 to 03-Jun-1996/23:59:59 B Jones Initial fit per unit weight = 0.442085 Jy, sum of weights = 7.1308e+06 B Jones Final fit per unit weight = 0.0154699 Jy after 7 iterations Storing G matrix in table n4038D.gcal Caught exception: Table n4038D.gcal cannot be created; it is in use in another process Caught an exception! Event type=run exception=(/aips++/weekly/code/trial/implement/MeasurementComponents/DOcalibrater.cc : 401) Failed AlwaysAssertExit retval LOG> messages in glish window: - cal.solve(); : Caught an exception! Event type=run exception=(/aips++/weekly/code/trial/implement/MeasurementComponents/DOcalibrater.cc : 401) Failed AlwaysAssertExit retval File: servers.g, Line 931 Stack: .(), calibrater.g line 264 .() BUG> AIPS++ defect AOCso03051: tries to re-write G matrix - cal.done(); JEH> OK, enough for today. Need some resolution on the new bugs before JEH> I can proceed. JEH> 10/03/01 Going to repeat, but opening and closing the cal tool as JEH> in the script. - include 'calibrater.g'; - cal:=calibrater('n4038D.ms'); - cal.setdata(msselect='FIELD_ID<=2'); - cal.setsolve(type='G',t=900.0,preavg=0.0,refant=24,table='n4038D.gcal'); - cal.solve(); - cal.plotcal(tablename='n4038D.gcal') - cal.fluxscale(tablein='n4038D.gcal', tableout='n4038D.gcalout', reference='1331+305', transfer=['1130-148']); - cal.plotcal(tablename='n4038D.gcalout') - cal.done(); - cal:=calibrater('n4038D.ms'); - cal.setdata(msselect='FIELD_ID==1'); - cal.setapply(type='G',t=0.0,table='n4038D.gcalout'); - cal.setsolve(type='B',t=86400.0,preavg=86400.0,refant=24,table='n4038D.bcal'); - cal.solve(); LOG> Starting calibrater::solve LOG> Initializing solvable bandpass terms (B-matrix) LOG> For interval of 86400 seconds, found 1 slots LOG> Loading G table from n4038D.gcalout LOG> Solving for B LOG> B Jones Slot=1, 1331+305, spw=1: 03-Jun-1996/24:00:00 to 03-Jun-1996/23:59:59 LOG> B Jones Initial fit per unit weight = 1.49751 Jy, sum of weights = 7.1308e+06 LOG> B Jones Final fit per unit weight = 0.0139433 Jy after 6 iterations LOG> Storing B matrix in table n4038D.bcal LOG> Finished calibrater::solve LOG> 16.21 real 13.51 user 1.59 system - cal.plotcal(tablename='n4038D.bcal'); - cal.done(); JEH> Well, I'll be damned. It worked. Lets go ahead and repeat it JEH> withOUT the cal.done() in the middle and verify that I get JEH> the same error: Yep, exact same. So that's stupid. Updated my JEH> bug with this info. JEH> Georges' email said to restrict the G solution to the good JEH> channels. I'll try that to see if I can get the getjy fluxes JEH> to work. First, delete the old tables. - include 'calibrater.g'; - cal:=calibrater('n4038D.ms'); #use a psuedo-continuum from channels 11-117 for calculating the G Jones - cal.setdata(mode='channel', nchan=1, start=11, step=96, msselect='FIELD_ID<=2'); - cal.setsolve(type='G',t=900.0,preavg=0.0,refant=24,table='n4038D.gcal'); - cal.solve(); - cal.plotcal(tablename='n4038D.gcal') - cal.fluxscale(tablein='n4038D.gcal', tableout='n4038D.gcalout', reference='1331+305', transfer=['1130-148']); #LOG> Set the flux density scale using reference calibrators #LOG> Flux density for 1130-148 (spw=1) is: 2.171 +/- 0.002 Jy #JEH> So that didn't do any good. - cal.plotcal(tablename='n4038D.gcalout') - cal.done(); #now calculate bandpass, after first applying the G Jones - cal:=calibrater('n4038D.ms'); - cal.setdata(mode='channel', nchan=127, start=1, step=1, msselect='FIELD_ID==1'); - cal.setapply(type='G',t=0.0,table='n4038D.gcalout'); - cal.setsolve(type='B',t=86400.0,preavg=86400.0,refant=24,table='n4038D.bcal'); - cal.solve(); #JEH> exact same result as before. - cal.plotcal(tablename='n4038D.bcal'); #JEH> Lets look at more than this: - cal.plotcal(tablename='n4038D.bcal',multiplot='T',nx=3,ny=3); - cal.plotcal(tablename='n4038D.bcal',antennas=[1,2,3,4,5,6,7,8,9], multiplot='T',nx=3,ny=3); - cal.plotcal(tablename='n4038D.bcal',antennas=[9], multiplot='T',nx=3,ny=3); : Nothing to plot for specified selection. File: note.g, Line 58 Stack: throw(), calutil.g line 796 .(), calibrater.g line 322 .() - cal.done(); #JEH> Couldnt get this to work. Tried many antenna>1; same error. #BUG> AIPS++ defect AOCso03052: only appears to be one antennae in bpcal #JEH> More tomorrow. #JEH> No replies - everyone at ADASS. So I'll just continue on with #JEH> the script. This is not how I would have proceeded - I would #JEH> have first mapped ch0 to make sure calibration was ok, and #JEH> after then mapped, but not cleaned. Then subtract a mean of the #JEH> dirty image from the dirty cube. But I'll proceed with this #JEH> just b/c I am bored and want to see something. - cal:=calibrater('n4038D.ms'); - cal.setdata(msselect='FIELD_ID IN [1:3]'); - cal.setapply(type='G',t=0.0,table='n4038D.gcalout',select='FIELD_ID==2'); #JEH> Why field_id==2??? - cal.setapply(type='B',t=0.0,table='n4038D.bcal'); - cal.correct(); - cal.done(); #LOG> Starting calibrater::correct #LOG> Loading B table from n4038D.bcal #LOG> Loading G table from n4038D.gcalout #LOG> Finished calibrater::correct #LOG> 47.35 real 35.1 user 3.82 system - include 'imager.g' - imagr:=imager('n4038D.ms'); - imagr.setdata(fieldid=3, mode='channel', nchan=120, start=1, step=1, spwid=1); #LOG> Starting imager::setdata #LOG> Selecting data #LOG> Performing selection on MeasurementSet #LOG> Selecting on field and spectral window ids #LOG> No selection string given #LOG> By selection 77213 rows are reduced to 51394 #LOG> Selecting 120 channels, starting at visibility channel 1 stepped by 1 #LOG> Finished imager::setdata #LOG> 0.73 real 0.24 user 0.15 system - imagr.setimage(nx=256, ny=256, cellx='10arcsec', celly='10arcsec', stokes='I', fieldid=3, spwid=1, start=1, step=1, nchan=120, mode='channel'); - imagr.weight(robust=1); #LOG> Starting imager::weight #LOG> Weighting MS: IMAGING_WEIGHT column will be changed #LOG> Uniform weighting: sidelobes will be suppressed over full image #LOG> : using 256 pixels in the uv plane #LOG> Sum of weights = 1384 #LOG> Finished imager::weight #LOG> 8.1 real 6.08 user 1.12 system - imagr.make('n4038D.model'); #LOG> Starting imager::make #LOG> Making empty image: n4038D.model #LOG> Image spectral coordinate: 120 channels, starting at visibility channel 1 stepped by 1 #LOG> Frequency = 1.41109, channel increment = 2.44141e-05GHz #LOG> Rest frequency is 1e-09GHz #LOG> Image polarization = Stokes I #LOG> Finished imager::make #LOG> 4 real 0.86 user 0.66 system - imagr.clean(algorithm='clark', niter=3000, gain=0.1, threshold='0.0Jy', model='n4038D.model', image='n4038D.restored', residual='n4038D.residual'); #LOG> Starting imager::clean #LOG> Cleaning images #LOG> Found 1 specified model images #LOG> Clean gain = 0.1, Niter = 3000, Threshold = 0 Jy, Algorithm = clark #LOG> Fourier transforms will use image centers as tangent points #LOG> made generic SkyEquation #LOG> Created Sky Equation #LOG> Starting deconvolution #LOG> Transforming I only #LOG> Processing channel 1 of 120 #LOG> Initial maximum residual: 0.154387 #LOG> Iteration: 23, Maximum residual=0.122566 #LOG> Iteration: 409, Maximum residual=0.0947001 #LOG> Iteration: 1681, Maximum residual=0.0692252 #LOG> Iteration: 3000, Maximum residual=0.0817687 #LOG> Clean used 3000 iterations to get to a max residual of 0.0817687 #LOG> Processing channel 2 of 120 #LOG> Initial maximum residual: 0.123616 #LOG> Iteration: 9, Maximum residual=0.0806369 #LOG> Iteration: 57, Maximum residual=0.0403004 #LOG> Iteration: 545, Maximum residual=0.0330554 #LOG> Iteration: 1432, Maximum residual=0.0169165 #LOG> Iteration: 3000, Maximum residual=0.0215636 #LOG> Clean used 3000 iterations to get to a max residual of 0.0215636 #LOG> Processing channel 3 of 120 #LOG> Initial maximum residual: 0.115991 #LOG> Iteration: 8, Maximum residual=0.073574 #LOG> Iteration: 54, Maximum residual=0.0370286 #LOG> Iteration: 289, Maximum residual=0.0191365 #LOG> Iteration: 1695, Maximum residual=0.0147843 #LOG> Iteration: 3000, Maximum residual=0.0133958 #LOG> Clean used 3000 iterations to get to a max residual of 0.0133958 #LOG> ... #LOG> Processing channel 119 of 120 #LOG> Initial maximum residual: 0.117546 #LOG> Iteration: 12, Maximum residual=0.0710804 #LOG> Iteration: 78, Maximum residual=0.030747 #LOG> Iteration: 802, Maximum residual=0.0218382 #LOG> Iteration: 2844, Maximum residual=0.010807 #LOG> Iteration: 3000, Maximum residual=0.00774578 #LOG> Clean used 3000 iterations to get to a max residual of 0.00774578 #LOG> Processing channel 120 of 120 #LOG> Initial maximum residual: 0.108535 #LOG> Iteration: 13, Maximum residual=0.0662591 #LOG> Iteration: 96, Maximum residual=0.0300094 #LOG> Iteration: 1041, Maximum residual=0.0195314 #LOG> Iteration: 3000, Maximum residual=0.0137447 #LOG> Clean used 3000 iterations to get to a max residual of 0.0137447 #LOG> Successfully deconvolved image #LOG> Fitted beam used in restoration: 38.4494 by 20.1039 (arcsec) at pa 76.0717 (deg) #LOG> Finished imager::clean #LOG> 407.87 real 292.78 user 39.52 system - imagr.restore(model='n4038D.model', image='n4038D.restored2'); #LOG> Starting imager::restore #LOG> Restoring 1 models #LOG> Using previous beam fit #LOG> made generic SkyEquation #LOG> Transforming I only #LOG> Transforming I only #LOG> Finished imager::restore #LOG> 202.77 real 159.13 user 18.04 system - imagr.done(); # # Subtract mean continuum channels (5-23, 105-115) using image tool # - include 'image.g'; # Mask bad channels (1-5, 115-120: - restim:=image('n4038D.restored2'); - restim.set(pixelmask=F,region=drm.box([1,1,1,1],[256,256,1,5])); - restim.set(pixelmask=F,region=drm.box([1,1,1,115],[256,256,1,120])); # Copy restored image for in-place continuum subtraction: - csubim:=imagefromimage(outfile='n4038D.final', infile='n4038D.restored2'); # Define region of image which is continuum (5-23, 105-115) - cregion1:=drm.box([1,1,1,5],[256,256,1,23]); - cont1:=restim.subimage(outfile='cont1.im',region=cregion1); - cregion2:=drm.box([1,1,1,105],[256,256,1,115]); - cont2:=restim.subimage(outfile='cont2.im',region=cregion2); - cont:=imageconcat(infiles="cont1.im cont2.im",relax=T); # Form integrated continuum as moment=-1 of continuum channels - cont.moments(outfile='n4038D.cont',moments=-1,axis=4); # Subtract continuum from cube: - csubim.calc('$csubim - n4038D.cont') # View final image - csubim.view(); #BUG> AIPS++ defect AOCso03053: get rid of all that black area! #BUG> AIPS++ defect AOCso03054: Search string doesn`t work in obvious way #JEH> Compare to AIPS cube: - imagefromfits('imaips', '/export/data_1/jhibbard/DATA/n4038/HIdata/N4038Dr+1.LMV.fits'); #LOG> Created image of shape[111, 126, 86, 1] #LOG> Copy FITS file to '/export/data_1/jhibbard/aips++/imaips' #LOG> All pixels fit in memory (1202796 pixels). # frames 1-5 & 115-120 are blank; this wasn't true of prior # example. # End #============================================================= #JEH> FINALLY back on the case 11/07/01 (HST Review panel intervened). #JEH> A number of bugs have been resolved. Most notably, msplot #JEH> works in raster mode, and george sent a script for plotting #JEH> the BP for all antennaes. So lets try both of those. #JEH> First, I removed the previous flag and all msplot directories. #JEH> Using weekly - include 'msplot.g' - mymsplot:=msplot('n4038D.ms', edit=T, flagfile="n4038Dflag01.flag"); #GUI> Show: Plot X vs Y, #GUI> Select: X axis=uvdist, Yaxis=data #GUI> Data selection: #GUI> Field Id=1 #GUI> Spectral selection: #GUI> List of spectral windows [1] #GUI> Number of channels 127 #GUI> Start channel 1 #GUI> Width of channel 1 #GUI> Channel increment 1 #GUI> #JEH? Don't forget to set the all channels, all correlations flags. #JEH> X=uvdist, Y=data #JEH> X=time, Y=data #JEH> X=time, Y=data, dataselect fieldid, iterate on antennae 1 #JEH> One last thing; lets try the tvflag type plot: #GUI> Show: display data as an image #GUI> Xaxis=time, Y-axis=ifr_number, Image pixels=corrected #GUI> Data selection: all unset #GUI> Spectral selection: #GUI> List of spectral windows [1] #GUI> Number of channels 1 #GUI> Start channel 11 #GUI> Width of channel 96 #GUI> Channel increment 1 #GUI> #JEH> Pretty groovy. But somethings up with my calibration. #JEH> Lets look at the Bpass table, using #JEH> georges new script for plotting bpass #JEH> (Modified slightly by me) include 'table.g' include 'pgplotter.g' t:=table('n4038D.bcal') ants:=t.getcol('ANTENNA1') + 1; Nants:=shape(ants) gain:=t.getcol('GAIN'); # the indicies on the gain array are # [pol1,pol2,spw,channel,antenna] # pol1 and pol2 are 1 or 2, and are # meaningful for B tables only if pol1=pol2 solok:=t.getcol('SOLUTION_OK'); # solok contains a 3D array of booleans # indicating which [spw,channel,antenna] # indicies are ok function plot1BP(i_ant,xylimits,gaintable,okarray) { Nchan:=shape(gaintable)[4] Npol:=shape(gaintable)[1] mypg.env(xylimits[1],xylimits[2], xylimits[3], xylimits[4], 0,0); if (any(okarray[1,,i_ant])) # does any channel for this antenna have ok solutions? { for (pol in 1:Npol) { x:=seq(Nchan) # list of channel numbers, from 1..Nchan y:=abs(gain[pol,pol,1,,i_ant]) # mask out the bad channels: x:=x[okarray[1,,i_ant]]; y:=y[okarray[1,,i_ant]]; mypg.sci((i_ant+pol-1)%10+2) # use colors 2-11 mypg.pt(x,y,17) # plot y vs. x with symbol 17 #make and plot antenna labels: xtxt:= (i_ant/Nants)*limits[2]; ytxt:= limits[4]*0.95; if (pol==2) ytxt:=limits[4]*0.9; mypg.ptxt(xtxt,ytxt,0,0,spaste(ants[i_ant])); } } } mypg:=pgplotter() mypg.ask:=T; #this doesn't seem to work - I want it to ask before #plotting the next bp. limits:= [0,Nchan+1,0,max(abs(gain))*1.1] for (i in 1:Nants-1) plot1BP(i,limits,gain,solok) limits := [2,125,0.5,0.7] for (i in 1:Nants-1) plot1BP(i,limits,gain,solok) # reset color to white: mypg.sci(1) #JEH> Mean BP is ~0.6, so something is screwy. #JEH> I'd better start again, and do parallel processing in aips. #JEH> End of this file. See n4038D.txt2 for continuation. #=============================================================