In [1]:
from readRinexObs import rinexobs
import matplotlib.pyplot as plt
%matplotlib inline
import numpy as np
from matplotlib.dates import YearLocator, MonthLocator, DateFormatter
from pandas import DataFrame, Panel4D, read_hdf
from glob import glob
from datetime import datetime
import time

In [3]:
files = glob("/home/greg/Documents/Summer Research/rinex files/mah*")

svn = 23

convertL1ToMeters = 3.0E8/(154.0*10.23E6)
convertL2ToMeters = 3.0E8/(120.0*10.23E6)
f2f1Factor = 1.545727
convertMetersToTEC = 6.158

start = datetime(year=2015,month=10,day=7,hour=6,minute=0,second=49)
end = datetime(year=2015,month=10,day=7,hour=7,minute=7,second=30)

dat = {}

t=time.time()
for file in files:  
    d = rinexobs(file)
    dat[file.split('/')[-1]] = d[start:end,svn,:,:]
data = Panel4D(dat)
print("Seconds to read all files {0:.2f}".format(time.time()-t))
data.to_hdf('svn23.h5','data',format='table')


/home/greg/Documents/Summer Research/rinex files/mah62800.15o is a RINEX 2.11 file, 31649.899 kB.
30.61 seconds for _block2df
17.10 seconds for panel assignments
finished in 50.99 seconds
/home/greg/Documents/Summer Research/rinex files/mah32800.15o is a RINEX 2.11 file, 38072.718 kB.
31.76 seconds for _block2df
17.07 seconds for panel assignments
finished in 51.97 seconds
/home/greg/Documents/Summer Research/rinex files/mah52800.15o is a RINEX 2.11 file, 35699.333 kB.
32.45 seconds for _block2df
16.99 seconds for panel assignments
finished in 52.57 seconds
/home/greg/Documents/Summer Research/rinex files/mah72800.15o is a RINEX 2.11 file, 35799.149 kB.
31.43 seconds for _block2df
17.06 seconds for panel assignments
finished in 51.96 seconds
/home/greg/Documents/Summer Research/rinex files/mah42800.15o is a RINEX 2.11 file, 37584.842 kB.
35.67 seconds for _block2df
20.14 seconds for panel assignments
finished in 59.28 seconds
/home/greg/Documents/Summer Research/rinex files/mah92800.15o is a RINEX 2.11 file, 36138.917 kB.
35.07 seconds for _block2df
20.08 seconds for panel assignments
finished in 58.55 seconds
/home/greg/Documents/Summer Research/rinex files/mah22800.15o is a RINEX 2.11 file, 15777.532 kB.
20.55 seconds for _block2df
11.86 seconds for panel assignments
finished in 34.67 seconds
/home/greg/Documents/Summer Research/rinex files/mah82800.15o is a RINEX 2.11 file, 38070.627 kB.
36.80 seconds for _block2df
22.67 seconds for panel assignments
finished in 64.86 seconds
Seconds to read all files 425.33
---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in _create_storer(self, group, format, value, append, **kwargs)
   1168             try:
-> 1169                 return globals()[_STORER_MAP[pt]](self, group, **kwargs)
   1170             except:

KeyError: 'ndim'

During handling of the above exception, another exception occurred:

TypeError                                 Traceback (most recent call last)
<ipython-input-3-64a056fd2d17> in <module>()
     19 data = Panel4D(dat)
     20 print("Seconds to read all files {0:.2f}".format(time.time()-t))
---> 21 data.to_hdf('data.h5','data')

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/core/generic.py in to_hdf(self, path_or_buf, key, **kwargs)
   1094 
   1095         from pandas.io import pytables
-> 1096         return pytables.to_hdf(path_or_buf, key, self, **kwargs)
   1097 
   1098     def to_msgpack(self, path_or_buf=None, encoding='utf-8', **kwargs):

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in to_hdf(path_or_buf, key, value, mode, complevel, complib, append, **kwargs)
    258         with HDFStore(path_or_buf, mode=mode, complevel=complevel,
    259                       complib=complib) as store:
--> 260             f(store)
    261     else:
    262         f(path_or_buf)

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in <lambda>(store)
    253         f = lambda store: store.append(key, value, **kwargs)
    254     else:
--> 255         f = lambda store: store.put(key, value, **kwargs)
    256 
    257     if isinstance(path_or_buf, string_types):

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in put(self, key, value, format, append, **kwargs)
    824             format = get_option("io.hdf.default_format") or 'fixed'
    825         kwargs = self._validate_format(format, kwargs)
--> 826         self._write_to_group(key, value, append=append, **kwargs)
    827 
    828     def remove(self, key, where=None, start=None, stop=None):

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in _write_to_group(self, key, value, format, index, append, complib, encoding, **kwargs)
   1244 
   1245         s = self._create_storer(group, format, value, append=append,
-> 1246                                 encoding=encoding, **kwargs)
   1247         if append:
   1248             # raise if we are trying to append to a Fixed format,

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in _create_storer(self, group, format, value, append, **kwargs)
   1169                 return globals()[_STORER_MAP[pt]](self, group, **kwargs)
   1170             except:
-> 1171                 error('_STORER_MAP')
   1172 
   1173         # existing node (and must be a table)

/home/greg/anaconda3/lib/python3.5/site-packages/pandas/io/pytables.py in error(t)
   1134                 "cannot properly create the storer for: [%s] [group->%s,"
   1135                 "value->%s,format->%s,append->%s,kwargs->%s]"
-> 1136                 % (t, group, type(value), format, append, kwargs)
   1137             )
   1138 

TypeError: cannot properly create the storer for: [_STORER_MAP] [group->/data (Group) '',value-><class 'pandas.core.panelnd.Panel4D'>,format->fixed,append->False,kwargs->{'encoding': None}]

In [2]:
data = read_hdf('svn23.h5','data')
fig = plt.figure(figsize=(16,16))
ax1 = plt.subplot(211)
fmt = DateFormatter('%H:%M:%S')
ax1.xaxis.set_major_formatter(fmt)
ax1.autoscale_view()
plt.title('TEC: SVN 23 from different receivers')
plt.xlabel('Time')
plt.ylabel('TECu')


convertL1ToMeters = 3.0E8/(154.0*10.23E6)
convertL2ToMeters = 3.0E8/(120.0*10.23E6)
f2f1Factor = 1.545727
convertMetersToTEC = 6.158
handles=[]

for site in data.labels:
    
    lli = data[site,:,['L1','L2','C1','P2'],'lli']
    lli[np.isnan(lli)]=0
    lli=lli.astype(int)
    llimask = np.logical_or.reduce(lli%2)
        
    diffRange = data[site,:,'P2','data'] - data[site,:,'C1','data']
    diffRange = diffRange * convertMetersToTEC * f2f1Factor
    phase = f2f1Factor*(data[site,:,'L1','data']*convertL1ToMeters - data[site,:,'L2','data']*convertL2ToMeters)
    phase = phase * convertMetersToTEC
    
    avg = np.average((phase-diffRange)[np.logical_not(np.isnan(phase-diffRange))])
    TEC = phase - avg
    
    if(site!='mah92800.15o' and site!='mah72800.15o'): 
        handles.append(plt.plot(TEC[1100:1250],'.',label=site[:4]))
        plt.plot(TEC[1100:1250][llimask[1100:1250]],'kx',markersize=10,label='')
plt.legend()
plt.grid()
plt.show()