################### ANSWERS FOR THE INTERMEDIATE EXERCISE ######################
#import the packages we will need for this exercise 
import numpy as np
import matplotlib.pyplot as plt

#print the header for our table:
print 'SNR\tRMS(std)\tRMS(calc)' #the \t represents a tab in the string

#read in the file
f =open('intermediate.txt','r')

###################### LOOP OVER EACH LINE IN THE FILE #########################
#Note that each line contains one gaussian time series
for line in f.readlines():

################### CONVERT THE DATA TO A NUMPY ARRAY ##########################
	#strip the new line character (\n) off the end of the line, then split
	#the line into a list of entries
	data = line.strip().split(',')
	#convert the list into a numpy array (and convert all the entries from strings to floats)
	dataarray = np.array(data, dtype='float')
	
################################ MAKE A PLOT ###################################
	#determine the number of entries in your array, and then plot your data,
	#one data point per entry.
	Nentries = len(dataarray) 
	plt.plot(np.arange(0,Nentries,1),dataarray,'.') 
	#Note: arange generates a numpy array with syntax (first entry,length,step).
	#Further note: the '.' makes the it plot points - you can also use plt.scatter()
	#To display the plot to the screen we would use plt.show() here. However,
	#We're going to plot all the data on top of each other and show at the end,
	#so we don't use plt.show() right now 

####################### CALCULATE SIGNAL TO NOISE RATIO ########################	
	#find the peak flux
	peakflux = np.max(dataarray)

	#here are two different ways of calculating the RMS.
	#Note that we cheat here: we know the pulse is in the first half of the
	#time series, so we just calculate the rms noise on the second half
	rmsflux = np.std(dataarray[Nentries/2:])
	rmsflux2 = np.sqrt(np.mean(dataarray[Nentries/2:]**2))

	#print out your signal to noise calculation
	print '%(SN).2f\t%(peak).2f\t%(rms).2f\t%(rms2).2f' %{'SN':peakflux/rmsflux, 'peak':peakflux, 'rms':rmsflux, 'rms2':rmsflux2}

################################### END LOOP ###################################

############################# DISPLAY THE PLOT #################################
plt.xlabel('time (s)')
plt.ylabel('Flux (mJy)')
plt.title("Michael's Face When:")
plt.show()