Import and Process SHG Data

In [ ]:

    

import numpy as np

import pandas as pd

import matplotlib.pyplot as plt

from scipy.optimize import curve_fit

import math


    

In [ ]:

    

#location of data files
inputpath = '/Users/geiger/Box Sync/Science/Data/SHG/2018/120918'

#path of output files
outputpath = '/Users/geiger/Box Sync/Science/Data/SHG/2018/120918/'


    

In [ ]:

    

#names of each file if have individual names
#names = ['hemi2pure1a','hemi2pure1b','hemi2salt1a','hemi2salt1b','hemi2pure2a','hemi2pure2b','hemi2salt2a','hemi2salt2b','hemi2pure3a','hemi2pure3b']

#names of each file if have sequential names
names = []
basename = 'pm'
for i in range(26):
    names.append(basename+str(i+1))


    

In [ ]:

    

#open file for writing scans to
fScans = open(outputpath+'scans.txt','w+')

#initialize data frames to hold data
countsA = pd.DataFrame()
countsB = pd.DataFrame()
pos = pd.DataFrame()

#go through each file
for name in names:
    #names of each file
    filename = inputpath + '/' + name + '.txt'
    #import countsA (signal),countsB (dark counts), and pos (stage position)
    countsA[name] = pd.read_csv(filename,sep='\t')['countsA']
    countsB[name] = pd.read_csv(filename,sep='\t')['countsB']
    pos[name] = pd.read_csv(filename,sep='\t')['stage']
    
#define function to turn time wave into scan by taking average of intervals
def findAverage(series):
    #set number of points per position here
    interval = 20
    reshaped = np.reshape(series.values,(int(len(series.values)/interval),interval))
    return pd.Series(np.mean(reshaped,1))

#apply function to time data to get the scan data
aveCountsA = countsA.apply(findAverage,axis=0)
aveCountsB = countsB.apply(findAverage,axis=0)
pos = pos.apply(findAverage,axis=0)
del countsA,countsB
    
#correct for dark counts
counts = aveCountsA.sub(aveCountsB)
del aveCountsA,aveCountsB

#create data frame to hold each scan and the position vector
data = counts.copy()
data.insert(0,'pos',pos[names[0]])

#write to file    
fScans.write(data.to_csv(sep='\t', index=False, header=True))
fScans.close()

#plot individual
plt.figure()
plt.title('Individual Scan Data')
for column in pos.columns:
    plt.plot(pos[column],counts[column],'.')


    

In [ ]:

    

#define fit func, same as IGOR
def sinFunc(x, y0, A, f, phi):
    return y0 + A*np.sin(f*x + phi)

#x values from which to plot fit function
xvalues = np.linspace(0,99.7,1000)

#data frame to hold initial fitting parameters
initFitParams = pd.DataFrame(columns=['name','y0','y0error','A','Aerror','f','ferror','phi','phierror'])
initFitValues = pd.DataFrame({'pos':xvalues})

#fit, going through column by column and storing fit parameters in initFits
for column in counts.columns:
    #calculate guesses for fit func
    y0guess = np.mean(counts[column])
    Aguess = (np.amax(counts[column])-np.amin(counts[column]))/2
    fguess = 0.05;
    phiguess = 0;
    guesses = [y0guess,Aguess,fguess,phiguess]
    
    #fit it
    popt, pcov = curve_fit(sinFunc,pos[column],
                           counts[column],p0=guesses)
    #calculate standard error
    pstd = np.sqrt(np.diag(pcov))
    
    #create row and append it to the dataframe
    tdf = pd.DataFrame({'name':[column],'y0':[popt[0]],'y0error':[pstd[0]],'A':[popt[1]],'Aerror':[pstd[1]],
                         'f':[popt[2]],'ferror':[pstd[2]],'phi':[popt[3]],'phierror':[pstd[3]]})
    initFitParams = initFitParams.append(tdf,ignore_index=True)
    
    #calculate fit and add it to fit values
    initFitValues[column] = sinFunc(xvalues,popt[0],popt[1],popt[2],popt[3])
    

#resort columns
columnTitles = ['name','y0','y0error','A','Aerror','f','ferror','phi','phierror']
initFitParams = initFitParams.reindex(columns=columnTitles)

#plot the initial fits
plt.figure()
plt.title('Init Fits')
for column in counts.columns:
    plt.plot(pos[column],counts[column],'.')
    plt.plot(xvalues,initFitValues[column])


    

In [ ]:

    

#write parameters to file
fInitFitParams = open(outputpath+'initFitParams.txt','w+')
fInitFitParams.write(initFitParams.to_csv(sep='\t', index=False, header=True))
fInitFitParams.close()

#write values to file
fInitFitValues = open(outputpath+'initFitValues.txt','w+')
fInitFitValues.write(initFitValues.round(3).to_csv(sep='\t',index = False, header=True))
fInitFitValues.close()


    

In [ ]:

    

#calculate average of f values, then period
fAve = initFitParams['f'].mean()
period = 2*np.pi/fAve

#calculate stdev
fStd = initFitParams['f'].std()
periodError = period*(fStd/fAve)

#print
print('f = '+'%.2f'%period+' +- '+'%.2f'%periodError)


    

In [ ]:

    

#data frame to hold final fitting parameters and values
finalFitParams = pd.DataFrame(columns=['name','y0','y0error','A','Aerror','f','ferror','phi','phierror'])
finalFitValues = pd.DataFrame({'pos':xvalues})

#go through each column
for column in counts.columns:
    #calculate guesses
    y0guess = np.mean(counts[column])
    Aguess = (np.amax(counts[column])-np.amin(counts[column]))/2
    phiguess = 0;
    guesses = [y0guess,Aguess,phiguess]

    #fit it, with f fixed
    popt, pcov = curve_fit(lambda x, y0, A, 
                           phi: sinFunc(x,y0, A, fAve, phi),
                           pos[column],counts[column],p0=guesses)
    
    #calculate standard error
    pstd = np.sqrt(np.diag(pcov))
    
    #create row and append it to the dataframe
    tdf = pd.DataFrame({'name':[column],'y0':[popt[0]],'y0error':[pstd[0]],'A':[popt[1]],'Aerror':[pstd[1]],
                         'f':[fAve],'ferror':[fStd],'phi':[popt[2]],'phierror':[pstd[2]]})                       
    finalFitParams = finalFitParams.append(tdf,ignore_index=True)
    
    #calculate fit and add it to fit values
    finalFitValues[column] = sinFunc(xvalues,popt[0],popt[1],fAve,popt[2])
             
#calculate phi in degrees
finalFitParams['phideg'] = np.degrees(finalFitParams['phi'])
finalFitParams['phidegerror'] = (finalFitParams['phierror']/finalFitParams['phi'])*finalFitParams['phideg']

#resort columns
columnTitles = ['name','y0','y0error','A','Aerror','f','ferror','phi','phierror','phideg','phidegerror']
finalFitParams = finalFitParams.reindex(columns=columnTitles)

#plot the final fits
plt.figure()
plt.title('Final Fits')
for column in counts.columns:
    plt.plot(pos[column],counts[column],'.')
    plt.plot(xvalues,finalFitValues[column])


    

In [ ]:

    

#write parameters to file
fFinalFitParams = open(outputpath+'finalFitParams.txt','w+')
fFinalFitParams.write(finalFitParams.to_csv(sep='\t', index=False, header=True))
fFinalFitParams.close()

#write values to file
fFinalFitValues = open(outputpath+'finalFitValues.txt','w+')
fFinalFitValues.write(finalFitValues.round(3).to_csv(sep='\t',index = False, header=True))
fFinalFitValues.close()


    

In [ ]:

    




    

In [ ]: