#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Mon Dec 16 11:40:15 2019

@author: dasanna
"""

import numpy as np
import tidynamics
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

plt.rcParams['figure.figsize'] = 6.4, 4.8
plt.rcParams['figure.subplot.bottom'] = 0.14
plt.rcParams['figure.subplot.left'] = 0.12



N = 3499
mean = np.zeros(N)
count = 0


for i in range(0,100):
    # Generate steps of value +/- 1
    filename="PARA_POSITIONS/pos"+str(i+1)+".dat"
    # Compute random walk position
    data = np.loadtxt(filename)
    msd =   tidynamics.msd(data)
    msd = msd[1:N//2]
    mean += tidynamics.msd(data)
    time=0.0065*np.arange(N)[1:N//2]
    count += 1


mean /= count
mean=mean[1:N//2]


def func(x, a):
    return a*x

popt, pcov = curve_fit(func,time,mean)

print(popt)

tarray = np.linspace(0.05, 15.0, num=100)

plt.xscale('log')
plt.yscale('log')

plt.xlim(0.05,15)
plt.ylim(0.1,10)

plt.xticks(fontsize=22)
plt.yticks(fontsize=22)

plt.xlabel(r'$t$ (s)',fontsize=22)
plt.ylabel(r'MSD ($\mu m^2$)',fontsize=20)

plt.plot(tarray,func(tarray,*popt),'k-',linewidth=2,label=r'$D\,t$')


plt.plot(time[::2],mean[::2],'ro',markersize=6)

plt.savefig('Fig2_c.pdf',bbox_inches='tight')

plt.show()
