This is a mess...

Nick nleioatt at gmail.com
Thu Jul 16 10:07:58 EDT 2009


I've been coding python for about a week now, and i'm trying to make
an object oriented version of a program i just wrote.  in this post is
the original program.  the next post will include the new programs and
the traceback.  i'm sure there are many underlying problems, but have
been stuck on this for 2 days now...

email me for input files since they are big and i don't see a way to
upload a file here(~5 mB)
this program works, and inputing the same file twice you should get a
"Q" values of ~1 (.99...)

##############CODE#####################

#9July09
#Compute the covariance overlap of two results
#Read in _U.asc files for both ensembles to find e-vectors
#Read in _ .asc files for both ensembles to find sqrt e-values
#strip out matrix info
#compute corvariance


import sys
import math
from math import sqrt


if len(sys.argv) != 3:
    print " "
    print "The corrent usage is 'python covariance.py file1 file2'"
    print "where the _U.asc and _s.asc will be appended when needed"
    print " "
    exit(1)


####Input the first prefix1_U.asc file:
####this is the first eigen vector matrix
prefix1 = sys.argv[1] + "_U.asc"
fileholder = open(prefix1)
text = fileholder.readlines()
fields = text[1].split()
num_rows = int(fields[1])
num_cols = int(fields[2])

U1_matrix = []
for line in text[2: num_rows+2]:
    fields = line.split()
    for i in range(len(fields)):
        fields[i] = float(fields[i])
    U1_matrix.append(fields)

####Input the second prefix2_U.asc file:
####this is the 2nd eigen vector matrix
prefix2 = sys.argv[2] + "_U.asc"
fileholder = open(prefix2)
text = fileholder.readlines()
fields = text[1].split()
num_rows = int(fields[1])
num_cols = int(fields[2])
#print "look here nick:", fields

U2_matrix = []
for line in text[2: num_rows+2]:
    fields = line.split()
    for i in range(len(fields)):
        fields[i] = float(fields[i])
    U2_matrix.append(fields)

####Then Read in the first eigen values
####1st 2 lines are header and need to be stripped away
prefix3 = sys.argv[1] + "_s.asc"
fileholder = open(prefix3)
text = fileholder.readlines()
fields = text[1].split()
num_vals1 = int(fields[1]) #add check to see if correct # of values
added
print "square if", len(U1_matrix), "=", len(U2_matrix), "=", len
(U1_matrix[1]), "=", len(U2_matrix[1])
#the list of sqrt e vals
sqrt_s1 = [float(line) for line in text[2: num_vals1+2]]
s1 = [float(line) * float(line) for line in text[2: num_vals1+2]]

####Finally, read in the 2nd set of eigen vals
prefix4 = sys.argv[2] + "_s.asc"
fileholder = open(prefix4)
text = fileholder.readlines()
fields = text[1].split()
num_vals2 = int(fields[1]) #add check to see if correct # of values
added
#the list of sqrt e vals
sqrt_s2 = [float(line) for line in text[2: num_vals1+2]]
s2 = [float(line) * float(line) for line in text[2: num_vals1+2]]


#=============================================================
####double summation (the 2nd term)
doublesum = 0.0
for i in range(len(U1_matrix[1])):
   #print "i = ", i
    v1 = U1_matrix[i]
    sum = 0
    for j in range(len(U2_matrix)):
         v2 = U2_matrix[j]
         dot = 0
         for k in range(len(U1_matrix)):
             dot += v1[k] * v2[k]
         root = sqrt_s1[i] * sqrt_s2[j]
         sum += root * dot * dot
    doublesum += sum


print "double sum: ", doublesum

####single summation (the 1st term and denominator)
singsum = 0.0
for k in range(len(s1)):
    singsum += s1[k] + s2[k]
print "sing sum:", singsum


####Put it all together
Q = 1 - sqrt(abs((singsum - (2.0*doublesum)) / singsum))
print "your Q:"
print Q

############################# end of code covariance.py




More information about the Python-list mailing list