In [32]:
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

In [93]:
def orders(data):
    S, T, L = 0.,0.,0.
    N = len(data)
    for d in data:
        S += np.cos(2.*(d[2]))
        T += np.sin(2.*(d[2]))
    S/=N
    T/=N
    L = np.sqrt(S*S + T*T)
    return L, S, T

In [106]:
path = "/home/michael/msc/mcmd/output/"
fnames = ["T", "X", "edgevar/edge15"]
snaps = range(0,4001,5)
L = np.zeros((len(fnames), len(snaps)))
S = np.zeros((len(fnames), len(snaps)))
T = np.zeros((len(fnames), len(snaps)))

In [107]:
print len(snaps)


801

In [108]:
iFile = 0
for f in fnames:
    fin = open(path+f, 'r')
    print "File", f
    iSnap = 0
    iRec = 0
    data = []
    record = 0 in snaps
    for line in fin.readlines():
        if line == "\n":
            if record:
                L[iFile][iRec], S[iFile][iRec], T[iFile][iRec] = orders(data)
                iRec+=1
            data = []
            if iSnap==snaps[-1]: break
            iSnap += 1
            record = iSnap in snaps
            if iSnap%(snaps[-1]/10) == 0: print iSnap
            continue
        elif record: 
            spt = line.split()
            spt = [float(x) for x in spt]
            data.append([spt[i] for i in [0,1,4,9,10,11,12,13,14,15,16]])

    fin.close()
    iFile += 1
print "Done"


File T
400
800
1200
1600
2000
2400
2800
3200
3600
4000
File X
400
800
1200
1600
2000
2400
2800
3200
3600
4000
File edgevar/edge15
400
800
1200
1600
2000
2400
2800
3200
3600
4000
Done

In [91]:
Ls,Ss,Ts = L,S,T

In [109]:
f, ax = plt.subplots(3,1)
for i in range(len(fnames)):
    im0 = ax[0].plot(snaps,L[i], label=fnames[i])
    im1 = ax[1].plot(snaps,S[i], label=fnames[i])
    im2 = ax[2].plot(snaps,T[i], label=fnames[i])

plt.legend(loc=7)
plt.gcf().set_size_inches(10,10)



In [110]:
f.savefig("TXD-compare.eps")

In [92]:
def orders_local(data):
    nbrthetas = []
    for d in data:
        select = []
        for dd in data:
            if dd[0] == d[0]: continue
            if dd[1] == d[1]:
                select.append(dd[2])
            for ni in d[3:]:
                if ni == -1.: continue
                if dd[1] == ni:
                    select.append(dd[2])
        nbrthetas.append([d[2], select])

    S, T, L = 0.,0.,0.
    N = 0
    for n in nbrthetas:
        for nn in n[1]:
            S += np.cos(2.*(n[0] - nn))
            T += np.sin(2.*(n[0] - nn))
            N+=1

    S/=N
    T/=N
    L = np.sqrt(S*S + T*T)
    return L, S, T