Commit 05b73ebc authored by Javier Quinteros's avatar Javier Quinteros
Browse files

Add input of network and channel code and save in I32 data type

MSEED record fo 4096
parent 0795fd9d
......@@ -279,7 +279,7 @@ class SDS(Archive):
daytrace = trace.slice(starttime=auxstart,
endtime=auxend,
nearest_sample=False)
daytrace.write(fout, format='MSEED', reclen=512)
daytrace.write(fout, format='MSEED', reclen=4096)
# Move to the next day
auxstart = copy(auxend)
......@@ -115,6 +115,12 @@ def main():
parser.add_argument('--decimate', type=int, choices=[1, 5],
help='Factor by which the sampling rate is lowered by decimation.',
default=1)
parser.add_argument('-N', '--network',
help='Network code to store in the miniseed header (default: "XX")',
default='XX')
parser.add_argument('-C', '--channel',
help='Channel code to store in the miniseed header (default: "FSF")',
default='FSF')
parser.add_argument('-o', '--outstruct', type=str, choices=dictarchive.keys(),
help=helparchive, default='StreamBased')
parser.add_argument('--metadata', action='store_true', default=False,
......@@ -139,6 +145,14 @@ def main():
logs.error('End time is smaller than Start time.')
return
if len(args.network) != 2:
logs.error('Network code must be two alphanumeric characters')
return
if len(args.channel) != 3:
logs.error('Channel code must be three alphanumeric characters')
return
td = TDMS(args.filename, directory=args.directory, iterate='M' if args.metadata else 'D',
chstart=args.chstart, chstop=args.chstop, chstep=args.chstep,
starttime=dtstart, endtime=dtend, decimate=args.decimate)
......
......@@ -26,6 +26,11 @@ from math import ceil
from numbers import Number
# Some examples with PDN_1km show that for 16 files of 44 MB (704 MB) we see the following
# With simple encoding of I16 and a record of 4096 bytes we use 720MB
# With STEIM2, a data type of I32 and a record of 4096 bytes we use 504 MB
# With STEIM2, a data type of I32, a decimation factor of 5, and a record of 4096 bytes we use 109 MB
def tup2time(fraction, seconds):
# logs = logging.getLogger('tup2time')
# logs.debug('seconds: %s' % seconds)
......@@ -68,7 +73,7 @@ class TDMS(object):
def __init__(self, filename, directory='.', chstart=0, chstop=None, chstep=1, channels=None,
starttime=None, endtime=None, iterate='D', decimate=1, firfilter='fir235',
loglevel='INFO'):
networkcode='XX', channelcode='FSF', loglevel='INFO'):
"""Initialize the TDMS object selecting the data, channels and decimation
:param filename: Experiment to read and process. Usually the first part of the filenames
......@@ -93,6 +98,10 @@ class TDMS(object):
:type decimate: int
:param firfilter: Filter to apply in case of decimation (fir235 is the only option)
:type firfilter: str
:param networkcode: Network code of the experiment. It has to be two alphanumeric characters
:type networkcode: str
:param channelcode: Channel code of the experiment. It has to be three alphanumeric characters
:type channelcode: str
:param loglevel: Verbosity in the output
:type loglevel: str
"""
......@@ -120,9 +129,20 @@ class TDMS(object):
logs.error('channels must be a list of numbers or None')
raise TypeError('channels must be a list of numbers or None')
if not isinstance(networkcode, str) or len(networkcode) != 2:
logs.error('Network code has to be two alphanumeric characters')
raise TypeError('Network code has to be two alphanumeric characters')
if not isinstance(channelcode, str) or len(channelcode) != 3:
logs.error('Channel code has to be three alphanumeric characters')
raise TypeError('Channel code has to be three alphanumeric characters')
self.__chstart = chstart
self.__chstop = chstop
self.__chstep = chstep
self.__networkcode = networkcode
self.__channelcode = channelcode
# channels has priority. Otherwise, chstart, chstop and chstep are used
if channels is not None and isinstance(channels, list):
......@@ -183,9 +203,6 @@ class TDMS(object):
# Keys will be channel number and the values will be np.arrays
self.__buffer = dict()
# Datatype of each channel
self.__datatypes = dict()
# Dictionary to save the metadata defined in the file
self.metadata = dict()
......@@ -230,6 +247,7 @@ class TDMS(object):
self.__endian = None
self.__datatype = None
self.__datatypesize = None
self.__outdatatype = None
self.numchannels = None
self.__samples = None
self.__samplestart = None
......@@ -433,8 +451,15 @@ class TDMS(object):
if datatype is None:
raise Exception('datatype definition not found in any channel!')
if self.__data2mask[datatype][0] == 'h':
if self.__data2mask[datatype][0] == 'b':
self.__datatype = '%ci1' % self.__endian
self.__outdatatype = '%ci4' % self.__endian
elif self.__data2mask[datatype][0] == 'h':
self.__datatype = '%ci2' % self.__endian
self.__outdatatype = '%ci4' % self.__endian
elif self.__data2mask[datatype][0] == 'i':
self.__datatype = '%ci4' % self.__endian
self.__outdatatype = '%ci4' % self.__endian
elif self.__data2mask[datatype][0] == 'f':
self.__datatype = '%cf4' % self.__endian
else:
......@@ -509,8 +534,8 @@ class TDMS(object):
# In the case that we iterate to get data
# If no decimation is needed
if self.__decimate == 1:
for info in self.__iter_data__():
yield info
for data, stats in self.__iter_data__():
yield (data.astype(self.__outdatatype), stats)
else:
# Use an input buffer to store the data coming in chunks from __iter_data__
inbuf = dict()
......@@ -572,6 +597,7 @@ class TDMS(object):
startofchunk = False
# Change the headers to reflect the decimation
# FIXME Decimation factor is hardcoded
outbuf[ch]['stats']['sampling_rate'] = stats['sampling_rate']/5
# print('outbuf new: %s' % outbuf[ch]['stats'])
......@@ -581,12 +607,14 @@ class TDMS(object):
# Convolution of inbuf with filter and then leave the last 234 values (len(filter)-1)
# This convolution is performed on a long array (inbuf) and a default filter with 235 components
nodecimation[ch] = np.convolve(inbuf[ch]['data'], self.__filter, 'valid').astype(data.dtype)
# We save the result of the convolution as int32 to be able to use STEIM2 later
nodecimation[ch] = np.convolve(inbuf[ch]['data'], self.__filter, 'valid').astype(self.__outdatatype)
logs.debug('filtered: %d components' % len(nodecimation[ch]))
logs.debug('filtered[%d][:11] %s' % (ch, nodecimation[ch][:11]))
logs.debug('filtered[%d][-11:] %s' % (ch, nodecimation[ch][-11:]))
# Check if we can copy as many components as a multiple of the decimation factor
# FIXME Decimation factor is hardcoded
leftover = len(nodecimation[ch]) % 5
logs.debug('filtered: leave %d components for next iteration %s' % (leftover, nodecimation[ch][-leftover:]))
if leftover:
......@@ -636,13 +664,14 @@ class TDMS(object):
data = self.__readdata(channels=self.__channels)
# Loop through channels
for ch in self.__channels:
stats = {'network': 'XX', 'station': '%05d' % ch, 'location': '',
'channel': 'FH1', 'npts': len(data[ch]),
stats = {'network': self.__networkcode, 'station': '%05d' % ch, 'location': '',
'channel': self.__channelcode, 'npts': len(data[ch]),
'sampling_rate': self.sampling_rate,
'starttime': UTCDateTime(self.__twstart),
'mseed': {'byteorder': self.__endian,
'dataquality': 'D',
'record_length': 512}}
'record_length': 4096,
'blkt1001': {'timing_quality': 100.0}}}
logs.debug('Data length: %d; First component: %s' % (len(data[ch]), data[ch][0]))
yield data[ch], stats
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment