#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Mon Jul  6 08:41:49 2020

@author: abramfarley
"""
#imports
from random import gauss
from random import seed
from pandas import Series
from pandas.plotting import autocorrelation_plot
import numpy as np
import matplotlib.pyplot as plt


#import timesynth as ts
#%%
# seed random number generator
seed(2)
# create white noise series
series = [gauss(0.0, 1.0) for i in range(10000)]
series = Series(series)
 #%%
print(min(series), max(series))

# summary stats
print(series.describe())
#%%

#need to 
series = series + abs(min(series))
print(series.describe())



#%%
series.plot()

#%%
plt.plot(series)


#%%
#make series orders of magnitude smaller to match K and domain 
# 10^0 -> 10^-5
series = series/20000
print(series.describe())

#%%
#########
#Harmonic signal
samplingFrequency   = 50  #50
samplingInterval       = 1 / samplingFrequency
beginTime           = 0
endTime = 200

signal1Frequency     = 0.5
signal2Frequency     = 1
signal3Frequency     = 5

time        = np.arange(beginTime, endTime, samplingInterval)

amplitude1 = 2*np.sin(2*np.pi*signal1Frequency*time)

amplitude2 = 2*np.sin(2*np.pi*signal2Frequency*time)

amplitude3 = 2*np.sin(2*np.pi*signal3Frequency*time)

figure, axis = plt.subplots(3, 1)

plt.subplots_adjust(hspace=1)

### signal 1
axis[0].set_title('Sine wave with a frequency of 1/86400 Hertz')

axis[0].plot(amplitude1)

axis[0].set_xlabel('Time (day)')

axis[0].set_ylabel('Amplitude')

axis[0].set_xlim(0,1000)


### signal 2
axis[1].set_title('Sine wave with a frequency of 0.1/86400 Hertz')

axis[1].plot(amplitude2)

axis[1].set_xlabel('Time')

axis[1].set_ylabel('Amplitude')

axis[1].set_xlim(0,1000)


### signal 3
axis[2].set_title('Sine wave with a frequency of 0.25/86400 Hertz')

axis[2].plot(amplitude3)

axis[2].set_xlabel('Time')

axis[2].set_ylabel('Amplitude')

axis[2].set_xlim(0,1000)




#%%
amplitude = amplitude1 + amplitude2 + amplitude3

figure, axis = plt.subplots(2, 1)

plt.subplots_adjust(hspace=1)

axis[0].set_title('Composition of Sine Waves')

axis[0].plot(time, amplitude)

axis[0].set_xlabel('Timestep')

axis[0].set_ylabel('Amplitude')

axis[0].set_xlim(0,20)


fourierTransform = np.fft.fft(amplitude)/len(amplitude)
fourierTransform = fourierTransform[range(int(len(amplitude)/2))]

tpCount     = len(amplitude)

values      = np.arange(int(tpCount/2))

timePeriod  = tpCount/samplingFrequency

frequencies = values/timePeriod

axis[1].set_title('Fourier transform depicting the frequency components')

 

axis[1].plot(frequencies, abs(fourierTransform))

axis[1].set_xlim(0,5.5)

axis[1].set_xlabel('Frequency')

axis[1].set_ylabel('Amplitude')

print((max(abs(fourierTransform))))
print(len(fourierTransform))

#%%
shifted_amplitude = amplitude + abs(min(amplitude)*2)
damped_amplitude = shifted_amplitude/1000   #4000

#%%
harmonic_series = Series(damped_amplitude)
print(harmonic_series.describe())




#%%
plt.plot(harmonic_series)
plt.xlim(0,1500)


#%%
#white noise to match harmonics' std and mean
mu, sigma = 0.010712, 0.002450
s = np.random.normal(mu, sigma, 10000)

noise_series = Series(s)
print(noise_series.describe())

#%%
plt.plot(noise_series, label = "noise")
plt.plot(harmonic_series, label = 'harmonic')
plt.title("Time Series Input")
plt.xlabel("Timestep")
plt.ylabel("Magnitude")
plt.xlim(0,1000)
plt.legend()
#%%
#save as csv file
np.savetxt("/Users/abramfarley/ParF/uahpc/matching_inputs/over0.4freq_hn/harmonic_input_july.csv", harmonic_series, delimiter=",")
np.savetxt("/Users/abramfarley/ParF/uahpc/matching_inputs/over0.4freq_hn/noise_input_july.csv", noise_series, delimiter=",")

#%%

 
#%%