In [1]:
import pandas
import seaborn as sns
from datetime import datetime
import matplotlib.pyplot as plt
In [2]:
%matplotlib inline
pandas.set_option("display.max_rows", 10)
sns.set(style="darkgrid")
In [3]:
file = '/Users/schriste/Desktop/AC_H1_EPM_201207.txt'
In [4]:
file
Out[4]:
In [5]:
ace_col_names = ['FP6P_.761-1.22MEV_IONS', 'FP7P_1.22-4.97MEV_IONS', 'UNC_FP6P_.761-1.22MEV_IONS', 'UNC_FP7P_1.22-4.97MEV_IONS']
names = ['date', 'hour'] + ace_col_names
ace_data = pandas.read_csv(file, skiprows=74, delim_whitespace=True, names = names)
In [6]:
ace_data = ace_data.truncate(after=len(ace_data)-5)
In [7]:
times = [datetime.strptime(t[0:-4], '%d-%m-%Y %H:%M:%S') for t in ace_data['date'] + ' ' + ace_data['hour']]
In [8]:
ace_data['times'] = times
ace_data = ace_data.set_index('times')
In [9]:
ace_data = ace_data.drop('hour',1)
ace_data = ace_data.drop('date',1)
In [10]:
for col in ace_col_names:
ace_data[col] = ace_data[col].convert_objects(convert_numeric=True)
In [11]:
ace_data.dtypes
Out[11]:
In [12]:
for col in ace_col_names:
ace_data[col][ace_data[col] < 0] = np.nan
In [13]:
palette = sns.color_palette()
In [14]:
plt.figure()
ace_data[ace_col_names[0]].plot()
ace_data[ace_col_names[1]].plot()
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)
plt.title("EPAM/ACE Electron Proton Alpha Monitor")
plt.ylabel("Ion Flux")
plt.show()
In [15]:
def read_iswa_enlil_dat(file, col_names):
names = ['year', 'month', 'day', 'hour', 'minute'] + col_names
data = pandas.read_csv(file, skiprows=2, names=names, delim_whitespace=True)
times_str = [str(int(t[1]['year'])) + '-' + str(int(t[1]['month'])).zfill(2) + '-' + str(int(t[1]['day'])).zfill(2) + ' ' + str(int(t[1]['hour'])).zfill(2) + ':' + str(int(t[1]['minute'])).zfill(2) for t in data.iterrows()]
times = [datetime.strptime(t, '%Y-%m-%d %H:%M') for t in times_str]
data['times'] = times
data = data.set_index('times')
data = data.drop('year',1)
data = data.drop('month',1)
data = data.drop('day',1)
data = data.drop('hour',1)
data = data.drop('minute',1)
return data
In [16]:
def read_iswa_data_dat(file):
data = pandas.read_csv(file, delim_whitespace=True)
times_str = [str(t[0]) + ' ' + str(t[1]['Timestamp'])[0:-2] for t in data.iterrows()]
times = [datetime.strptime(str_time, '%Y-%m-%d %H:%M:%S') for str_time in times_str]
data['times'] = times
data = data.set_index('times')
data = data.drop('Timestamp',1)
return data
In [17]:
file = '/Users/schriste/Desktop/20120712_193500_ENLIL_time_line.dat.txt'
In [18]:
col_names = ['B_enl','V_enl','n_enl','T_enl']
enlil_data = read_iswa_enlil_dat(file, col_names=col_names)
In [19]:
enlil_data
Out[19]:
In [20]:
plt.figure()
enlil_data.plot(subplots=True, figsize=(10, 10))
plt.show()
Deriving the Kp index from ENLIL data
In [21]:
kp90exp=9.5-np.exp(2.17676-(0.000052001)*((enlil_data.V_enl)**1.3333)*((enlil_data.B_enl)**0.6666)*(np.sin((np.pi/2)/2)**(8/3.)))
kp135exp=9.5-np.exp(2.17676-(0.000052001)*((enlil_data.V_enl)**1.3333)*((enlil_data.B_enl)**0.6666)*(np.sin((3*np.pi/4)/2)**(8/3.)))
kp180exp=9.5-np.exp(2.17676-(0.000052001)*((enlil_data.V_enl)**1.3333)*((enlil_data.B_enl)**0.6666)*(np.sin((np.pi)/2)**(8/3.)))
In [22]:
plt.figure()
kp90exp.plot()
kp135exp.plot()
kp180exp.plot()
plt.title("ENLIL Predicted KP Index")
plt.ylabel("KP")
plt.show()
In [23]:
file = '/Users/schriste/Desktop/20120712_193500_ENLIL_Kp_timeline.dat.txt'
In [24]:
kp_col_names = ['Kp_90', 'Kp_135', 'Kp_180']
enlil_kp_data = read_iswa_enlil_dat(file, col_names=kp_col_names)
In [25]:
enlil_kp_data
Out[25]:
In [26]:
plt.figure()
for col in kp_col_names:
enlil_kp_data[col].plot()
plt.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)
plt.title("ENLIL Predicted KP Index")
plt.ylabel("KP")
plt.show()
a current bug in our GOES data fetcher is that it cannot handle multiple days in the time range
In [27]:
from sunpy.lightcurve import GOESLightCurve
from sunpy.time import TimeRange
In [28]:
tr = TimeRange(ace_data.index[0].to_datetime(), ace_data.index[-1].to_datetime())
In [29]:
tr
Out[29]:
In [30]:
goes = GOESLightCurve.create(tr)
In [31]:
plt.figure()
ax = goes.data['xrsa'].plot()
plt.yscale('log')
ax = goes.data['xrsb'].plot()
plt.show()
note to self: had to manually delete blank lines at the end of file otherwise times have nan's and fails to parse
In [32]:
file = '/Users/schriste/Desktop/ISWA_GOES.txt'
In [33]:
iswa_goes = read_iswa_data_dat(file)
In [34]:
iswa_goes
Out[34]:
In [35]:
iswa_goes['Long_Wave'][iswa_goes['Long_Wave'] < 0] = np.nan
In [36]:
plt.figure()
iswa_goes.plot()
plt.yscale('log')
plt.title('GOES 1-8 A (ISWA)')
plt.ylabel('Watts')
plt.show()
In [37]:
file = '/Users/schriste/Desktop/ISWA_KP.txt'
In [38]:
iswa_kp = read_iswa_data_dat(file)
In [39]:
plt.figure()
iswa_kp.plot()
plt.title('KP (ISWA)')
plt.ylabel('KP Index')
plt.show()
In [40]:
iswa_kp
Out[40]:
In [41]:
plt.figure()
iswa_kp.plot()
kp135exp.plot()
kp180exp.plot()
kp90exp.plot()
plt.legend()
plt.show()
Documentation for the following file http://www.srl.caltech.edu/STEREO/Public/HET_public.html
In [116]:
file = '/Users/schriste/Desktop/BeH12Jul.1m.txt'
In [140]:
def read_stereob_impact(file):
names = ['0', 'year', 'month', 'day', 'hhmm',
'eflux_0.7-1.4', 'eflux_0.7-1.4_uncert',
'eflux_1.4-2.8', 'eflux_1.4-2.8_uncert',
'eflux_2.8-4.0', 'eflux_2.8-4.0_uncert',
'pflux_13.6-15.1', 'pflux_13.6-15.1_uncert',
'pflux_14.9-17.1', 'pflux_14.9-17.1_uncert',
'pflux_17.0-19.3', 'pflux_17.0-19.3_uncert',
'pflux_20.8-23.8', 'pflux_20.8-23.8_uncert',
'pflux_23.8-26.4', 'pflux_23.8-26.4_uncert',
'pflux_26.3-29.7', 'pflux_26.3-29.7_uncert',
'pflux_29.5-33.4', 'pflux_29.5-33.4_uncert',
'pflux_33.4-35.8', 'pflux_33.4-35.8_uncert',
'pflux_35.5-40.5', 'pflux_35.5-40.5_uncert',
'pflux_40.0-60.0', 'pflux_40.0-60.0_uncert',
'pflux_60.0-100.0', 'pflux_60.0-100.0_uncert']
data = pandas.read_csv(file, delim_whitespace=True, skiprows=24, dtype={'hhmm': np.dtype('str')}, header=None, names=names[0:40])
dates_str = [str(row[1][1]) + '-' + str(row[1][2]).zfill(1) + '-' + str(row[1][3]).zfill(2) + ' ' + row[1][4][0:2] + ':' + row[1][4][2:4] for row in data.iterrows()]
times = [datetime.strptime(str_time, '%Y-%b-%d %H:%M') for str_time in dates_str]
data = data.drop('year', 1)
data = data.drop('month', 1)
data = data.drop('hhmm', 1)
data = data.drop('day', 1)
data = data.drop('0', 1)
data['times'] = times
data = data.set_index('times')
return data
In [141]:
data = read_stereob_impact(file)
In [145]:
cols = data.columns
In [143]:
data
Out[143]:
In [150]:
e_cols = ['eflux_0.7-1.4', 'eflux_1.4-2.8', 'eflux_2.8-4.0']
p_cols = ['pflux_13.6-15.1', 'pflux_14.9-17.1', 'pflux_17.0-19.3', 'pflux_20.8-23.8', 'pflux_23.8-26.4',
'pflux_26.3-29.7', 'pflux_29.5-33.4', 'pflux_33.4-35.8', 'pflux_35.5-40.5', 'pflux_40.0-60.0',
'pflux_60.0-100.0', ]
In [151]:
plt.figure()
for col in e_cols:
data[col].plot()
plt.legend()
plt.title('Electrons')
plt.show()
plt.figure()
for col in p_cols:
data[col].plot()
plt.title('Protons')
plt.legend()
plt.show()
In [159]:
data.resample('T')
plt.figure()
for col in e_cols:
data[col].resample('H', how='mean').plot()
plt.legend()
plt.title('Electrons')
plt.show()
In [ ]: