source:
branches/hpc33/tutorials/tutorial_3.rst@
2431
Last change on this file since 2431 was 1636, checked in by , 15 years ago | |
---|---|
File size: 3.9 KB |
Tutorial 3 - Spectral-line gaussian component fitting and rms noise levels
.. sectionauthor:: Shari Breen
Files
- mon.asap: Data file
Data Log
- ON-OFF Position switching mode with Hobart
- 4 scans (2 OFF and 2 ON)
- Dual pol, 1 IF
Instructions
- Work through the list of commands given in the text file to calibrate the data and fit gaussian components to the emission. The commands should be typed into ASAP line-by-line until you get to the ADVANCED FITTING step which will require you to experiment with the parameters.
- ADDITIONAL exercise: Determine the rms noise in the spectrum using the create mask function you learned in the first section and the stat command specified in the text file.
Examples - Fitting gaussians
If you only see a file called mon.asap.tar, untar it before getting started
# Load the scantable for the source into ASAP mon = scantable("mon.asap") mon.summary() # Select the on-source scans source=mon.get_scan([1,3]) # Select the reference scans ref=mon.get_scan([0,2]) # Make the quotient spectra and plot it quot=quotient(source,ref) plotter.plot(quot) # Remove the baseline and plot quot.auto_poly_baseline() plotter.plot(quot) # Average the 2 scans together av=quot.average_time(align=True) plotter.plot(av) # Check the the rest frequency is set correctly and define unit as km/s print av.get_restfreqs() av.set_unit("km/s") plotter.plot(av) # Scale the two polarizations separately sel=selector() sel.set_polarisations(0) av.set_selection(sel) av.scale(50.6) sel.set_polarisations(1) av.set_selection(sel) av.scale(32.5) # unset the selection av.set_selection() # Average the polarizations together avpol=av.average_pol() # Plot the data plotter.plot(avpol) # Create a mask that contains the emission msk=avpol.create_mask([-10,30]) # Set up some fitting parameters f=fitter() f.set_scan(avpol,msk) f.set_function(gauss=2) # Fit a gaussian to the emission, plot the emission and fit f.fit() f.plot() # Save the plot plotter.save("monr2.eps") # Get the fit parameters f.get_parameters() # Plot the fit residuals f.plot(residual=True) # Plot each of the fitted gaussians separately, with or without the fit parameters overlaid f.plot(components=0,plotparms=True) f.plot(components=1)
ADVANCED FITTING
Sometimes it is necessary to force some of the fitting pa- rameters such as peak or FWHM. For example, if you want to fix the peak of the first gaussian component to 39 Jy you can use the following command:
f.set_gauss_parameters(39, 10.2, 1, peakfixed=1, component=0)
You can now experiment with forcing the values of the different parameters. To look at the help file on this function type help fitter.set_gauss parameters Additional Exercise: Determine the rms noise in the spectrum.
- Create a mask that excludes the emission (if you leave it in the rms noise will be significantly inflated), using the same command as before except this time you will have to specify two ranges of values either side of the emission.
- Find the rms noise of the spectrum:
avpol.stats(stat="rms",mask=mymask)
The stats function can be used to extract many other statistics from the data such as the maximum and minimum values, the median value and many more. Use help to find out how to extract these statistics from the spectrum.