Page 2 - Data processing script

def reports(output):
    with open(output) as f:
        count = 1                  #this defines the number of tallies in the output file
        x = 101                    #this is the first surface number
        r = str('surface  %d' % x) #r is a variable used to help find the point where the results start.
        d = 0 #initializes depth variable
        c = 0 #intilizes the c variable used in the while loop
        p = 0 #initializes variable for calculating percent change
        #The following while loop skips each line until it comes across the line defined in r
        Flux = []
        SEM = []
        Depth = []
        Change = []
        while count < 10:
               while c != r:
                   line = f.readline()
                   c = line.strip()
               line = f.readline() #reads the energy and SEM
               numbers = line.strip().split()
               Flux.append(numbers[0])
               SEM.append(numbers[1])
               Depth.append(d)
               p = '%g'%(100*((float(Flux[0])-float(Flux[-1])) / float(Flux[0]))) #calculates percent change
               Change.append(p)
               d += 4
               count += 1
               x += 4
               r = str('surface  %d' %x)
        with open('Report.tex', 'w') as table:
            table.write(START)
            for T in zip(Depth,Flux,Change,SEM,):
                table.write('%d & %9.3e & %5.1f & %6.4f\\\\\n'%tuple(float(x) for x in T))
            table.write('\\hline\n')
            table.write(END)

        print('Report.tex created, now command   pdflatex Report  .')


START = r'''
\documentclass[12pt]{report} %

\begin{document} %

MCNPX Simulation For XXXXX %

\begin{tabular}{|lrrr|} %
\hline %
Depth&\multicolumn{1}{c}{Flux}&\multicolumn{1}{c}{Reduction}&\multicolumn{1}{c|}{SEM}\\ %
&\multicolumn{1}{c}{$\frac{n}{cm^2s}$}&\multicolumn{1}{c}{\%}&\multicolumn{1}{c|}{\ }\\ \hline %
'''

END = r'''
\end{tabular}
\end{document}
'''

>>> from rnaa import reports
>>> reports('naa')
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "rnaa.py", line 22, in reports
    energy.append(numbers[0])
IndexError: list index out of range
>>> 

 while b<4097:
            line = f.readline() #reads the energy and SEM
            numbers = line.strip().split()
            energy.append(numbers[0])
            counts.append(numbers[1])
            SEM.append(numbers[2])
            b+=1

def reports(output):
    with open(output) as f:
        r = str('cell  400') #r is a variable used to help find the point where the results start.
        c = 0 #intilizes the c variable used in the while loop
        #The following while loop skips each line until it comes across the line defined in r
        import numpy as np
        import matplotlib.pyplot as plt
       
        while c != r:
            line = f.readline()
            c = line.strip()
        
        line = f.readline()
        line = f.readline()
        b=0
        energy = []
        counts = []
        SEM = []
        while b<4097:
            line = f.readline() #reads the energy and SEM
            numbers = line.strip().split()
            energy.append(numbers[0])
            counts.append(numbers[1])
            SEM.append(numbers[2])
            b+=1
           
        table = open('Report.txt', 'w')
        table.write('================================================================ \n')
        table.write('MCNPX NAA Simulation For XXXXX                                   \n')
        table.write('================================================================ \n')
        table.write('================================================================ \n')
        table.write('Energy(MeV)          Counts          SEM             \n')
        table.write('================================================================ \n')
        for (E,C,S,) in zip(energy,counts,SEM,):    ####### zip is roughly a matrix transposition.
            table.write(str(E).ljust(8) + '     ' +
                        str(C).ljust(10) + '     ' +
                        str(S) + '\n')
        print('Report created')
        table.close()
        
        
        
        plt.rcParams['legend.loc'] = 'best'
        plt.figure(1)
        plt.xlim(1e-6, 3)
        plt.autoscale(enable=True, axis='y', tight=None)
        plt.plot(energy, counts, 'b')
        plt.xlabel('Energy (MeV)')
        plt.ylabel('Counts')
        plt.title('Gamma Spectrum')
        
        print('plots created')
        plt.show()