Python assignments

From AstroEdWiki
Revision as of 20:48, 22 March 2018 by WikiSysop (talk | contribs)
Jump to navigation Jump to search

This page contains assignments for our Research Methods - Programming with Python short course.


Homework 1: An assignment for very simple Python

Create a program called "gamma.py" that reads a value from the keyboard and computes and prints its Gamma function. Try 3/2 and compare it to the (square root of pi) divided by 2.

Run and see that the program returns the same values for these two quantities (within the floating point precision of your computer).


Homework 2: An assignment for solving problems with Python

Version 1 for optical physicists

The website http://www.filmetrics.com/refractive-index-database offers access to tables of the refractive index of many common materials used in optics. You can navigate to the website and pick ones of interest to download them to your computer. For materials that do not absorb, the tables offered will be for the real part of the index, usually called n, while for materials like metals that absorb light the site offers the absorption coefficient which is the imaginary part of the index as well.

Download several of the tables that appear to be of interest and combine them into a master file that has the name of the material and the tabular information associated with the name. You may format the file in whatever way you find easy, with the goal to have a single file of indices of refraction you can edit with a text editor. As a suggestion, it may help to organize the file in a way that it can also be read into a spreadsheet program such as Libre Office (Open Office), Google Docs, or Excel. Files in which entries are systematically delimited by spaces, tabs, or commas will do this, and the files that are downloaded from the website are "tab delimited".

Python can read and write csv files with a built-in module. Here's how the reader works:

 import sys
 import csv
 if len(sys.argv) == 1:
   print ("")
   print ("Usage: csv_to_txt.py infile.csv outfile.txt")
   print ("")
   sys.exit("Read a csv spreadsheet and output as text\n")
 elif len(sys.argv) == 3:
   infile = sys.argv[1]
   outfile = sys.argv[2]
 else:
   print ("")
   print ("Usage: csv_to_txt.py infile.csv outfile.txt")
   sys.exit("Read a csv spreadsheet and output as text\n")
 outfp = open(outfile, 'w')
 infp = open(infile, 'r')
 reader = csv.reader(infp)
 nrow = 1
 for newrow in reader:
   print newrow
   newline = ','.join(newrow)+'\n'       
   print (nrow, newline)
   outfp.write(newline)
   nrow = nrow +1
 
 outfp.close()
 infp.close()
 exit()

The csv module returns a list with a string in newrow, and the code

 newline = ','.join(newrow)+'\n'

takes elements from the list and joins them separated by commas to make something we can write out or print. If there is only one string in newrow, there will be no commas added.

The write operation is simple because it has its own built-in routine to write a row, as you can see with this piece of code:

 outfp = open(outfile, 'w')
 writer = csv.writer(outfp)
 writer.writerow( ('File Name', 'Directory', 'Target', 'Observation Date', 'Exposure Time', 'Filter', 'Telescope', 'Instrument', 'Created Date') )

In this example, the row is a tuple of strings that are entries for that row. Of course, the strings can be variables which are modified in the program, such as

 writer.writerow( (filename, dirname, imtarget, imdateobs, imexptime, imfilter, imtelescope, iminstrument, imcreatedate) )


The assignment is to write a program that will read your file, and then prompt for a material name and a wavelength. It should display the real (and imaginary, if available) parts of the index of refraction and then prompt again for another entry. Pressing Enter without entering information should cause the program to exit cleanly. Have at least two materials in the file, one of which has absorption. The more materials, the better of course!

Submit both your data file and your program to the class upload site when you have it working.


Version 2 for astronomers

On the class server we have downloaded a copy of the SKY2000 star catalog from the Vizier service. You will find it here, along with a pdf file that describes the catalog, and a short version of the first 100 lines from it for experimenting with your code. There is a similar directory with the HD catalog, one of the smaller catalogs of bright stars which is incorporated into SKY2000. The HD catalog is an astronomer's standard reference for bright stars.

These catalogs are in plain text, formatted by column, and has one star per line. Each line contains specific information about the star, including cross-referenced catalog numbers, celestial coordinates, magnitudes and spectral types. Download the data files and explanatory text or pdf to your computer (right click in your browser, and save).

Use the ngc_reader.example program we have discussed in class as a model, and write a program that reads and parses this file to find entries based on some aspect of the catalog. The challenge for this assignment is to write a program that (ideally on the command line, but if necessary using the input inside the program) finds any entries for the star with the HD number 128620. A complication for this catalog is that the HD entries are not the first data column, so you will have to parse them from inside each line. One way to do it is with something like this, after using readlines() to read and store the catalog line-by-line in hd_lines:

i = 0
for entries in hd_lines:
  if i == 0:
    hd_catalog = {hd_lines[i][35:43].rstrip() : hd_lines[i]}
  else:
    hd_catalog[hd_lines[i][35:43].rstrip()]] = hd_lines[i]
  i = i + 1 

However, in the columns that follow the HD number there is extra information to designate multiple stars with components that share the HD number, or uncertainties in the naming.

The program should print the right ascension, declination, magnitude, and spectral type for the star extracted from the catalog entry.

Upload your finished working program to the class homework upload site.


Supplemental Homework A: An assignment to use Fourier Transforms

Use the numpy_fft.example as a model to write a program that will find the Fourier transform of an oscillator with two simultaneous frequencies and damping constants. Try one at 6 Hz, damping in 2 seconds, and another at 6.6 Hz, damping in 0.5 seconds, both with the same amplitude.

Rather than the single cosine of the example, use a sum of two damped cosines as input to the FFT. The other parts of the code should be mostly unchanged. You'll see beating between the two frequencies, and the effect of the damping on the linewidth.

Upload your finished working program to the class homework upload site.


Supplemental Homework B: An assignment to image a star

Use the display_gaussian_noise.example as a model to write a program that will create a color image of the diffraction pattern of a star seen with a telescope.

Use a 1000x1000 image array with each pixel the equivalent a small angle on the sky -- an image of a star. In the absence of diffraction all the light would be in the pixel at the center of the image. With diffraction, the light is spread over the image in a diffraction pattern for a circular aperture given by

A(theta) = (2 J_1(ka sin(theta) / ka sin(theta) )^2

where J_1 is a Bessel function of order 1, k is 2 pi / lambda, a is the radius of the aperture, and theta is the angle to the optical axis. In other words, theta is equivalent the distance from the center of the image to the center of each pixel on a scale of radians.

In SciPy, you can compute a Bessel function J_1 with simply

import scipy.special as sp
y = sp.j1(x)

where "x" is the argument and y is the value you want. Since this is a SciPy function, x can be a NumPy array and jn will return y as an array of Bessel function values. Notice that as x goes to zero, J_1(x)/x goes to 0.5 .

For this exercise try 0.01 arcseconds on the sky for each pixel so the 1000x1000 image will cover 10 arcseconds. (One arcsecond is 1/3600 degree or 4.848x10^-6 radians.) Use a telescope diameter of 0.5 meters, or a radius of 0.25 meters. Calculate the image in each of three wavelengths to approximate the eye's response to various wavelengths lambda (red: 650 nm, green: 520 nm, and blue: 450 nm), and make color images showing the diffraction patterns that theoretically limit the image quality of a 0.5 meter telescope.

If you have time, try generating images of two stars slightly separated by adding the offset diffraction patterns. You can show that they appear as one star when too close, and are clearly two stars when they are farther apart. However, you'll also see that this resolution depends on the color so it should make for interesting images.