From 354077bea65daa23b1b60bf465ee5db5e0de995c 2017-10-31 18:44:22 From: José Chávez Date: 2017-10-31 18:44:22 Subject: [PATCH] Merge branch 'v2.3' of http://jro-dev.igp.gob.pe/rhodecode/schain into v2.3 --- diff --git a/schainpy/CHANGELOG.md b/schainpy/CHANGELOG.md index 5f8cfa0..4f4d159 100644 --- a/schainpy/CHANGELOG.md +++ b/schainpy/CHANGELOG.md @@ -1,6 +1,9 @@ ## CHANGELOG: ### 2.3 +* Added support for Madrigal formats (reading/writing). +* Added support for reading BLTR parameters (*.sswma). +* Added support for reading Julia format (*.dat). * Added high order function `MPProject` for multiprocessing scripts. * Added two new Processing Units `PublishData` and `ReceiverData` for receiving and sending dataOut through multiple ways (tcp, ipc, inproc). * Added a new graphics Processing Unit `PlotterReceiver`. It is decoupled from normal processing sequence with support for data generated by multiprocessing scripts. diff --git a/schainpy/model/graphics/jroplot_data.py b/schainpy/model/graphics/jroplot_data.py index 71a9ccd..dccdf06 100644 --- a/schainpy/model/graphics/jroplot_data.py +++ b/schainpy/model/graphics/jroplot_data.py @@ -50,14 +50,15 @@ class PlotData(Operation, Process): __missing = 1E30 __attrs__ = ['show', 'save', 'xmin', 'xmax', 'ymin', 'ymax', 'zmin', 'zmax', - 'zlimits', 'xlabel', 'ylabel', 'cb_label', 'title', 'titles', 'colorbar', - 'bgcolor', 'width', 'height', 'localtime', 'oneFigure', 'showprofile'] + 'zlimits', 'xlabel', 'ylabel', 'xaxis','cb_label', 'title', + 'colorbar', 'bgcolor', 'width', 'height', 'localtime', 'oneFigure', + 'showprofile', 'decimation'] def __init__(self, **kwargs): Operation.__init__(self, plot=True, **kwargs) Process.__init__(self) - self.contador = 0 + self.kwargs['code'] = self.CODE self.mp = False self.data = None @@ -87,14 +88,14 @@ class PlotData(Operation, Process): self.ymin = kwargs.get('ymin', None) self.ymax = kwargs.get('ymax', None) self.xlabel = kwargs.get('xlabel', None) - self.__MAXNUMY = kwargs.get('decimation', 300) + self.__MAXNUMY = kwargs.get('decimation', 200) self.showSNR = kwargs.get('showSNR', False) self.oneFigure = kwargs.get('oneFigure', True) self.width = kwargs.get('width', None) self.height = kwargs.get('height', None) self.colorbar = kwargs.get('colorbar', True) self.factors = kwargs.get('factors', [1, 1, 1, 1, 1, 1, 1, 1]) - self.titles = ['' for __ in range(16)] + self.titles = kwargs.get('titles', []) self.polar = False def __fmtTime(self, x, pos): @@ -359,7 +360,7 @@ class PlotData(Operation, Process): if self.xaxis is 'time': dt = self.getDateTime(self.max_time) xmax = (dt.replace(hour=int(self.xmax), minute=59, second=59) - - datetime.datetime(1970, 1, 1)).total_seconds() + datetime.datetime(1970, 1, 1) + datetime.timedelta(seconds=1)).total_seconds() if self.data.localtime: xmax += time.timezone else: @@ -369,9 +370,8 @@ class PlotData(Operation, Process): ymax = self.ymax if self.ymax else numpy.nanmax(self.y) Y = numpy.array([10, 20, 50, 100, 200, 500, 1000, 2000]) - i = 1 if numpy.where(ymax < Y)[ - 0][0] < 0 else numpy.where(ymax < Y)[0][0] - ystep = Y[i - 1] / 5 + i = 1 if numpy.where(ymax-ymin < Y)[0][0] < 0 else numpy.where(ymax-ymin < Y)[0][0] + ystep = Y[i] / 5 for n, ax in enumerate(self.axes): if ax.firsttime: @@ -429,15 +429,15 @@ class PlotData(Operation, Process): if self.nrows == 0 or self.nplots == 0: log.warning('No data', self.name) fig.text(0.5, 0.5, 'No Data', fontsize='large', ha='center') + fig.canvas.manager.set_window_title(self.CODE) continue fig.tight_layout() fig.canvas.manager.set_window_title('{} - {}'.format(self.title, self.getDateTime(self.max_time).strftime('%Y/%m/%d'))) - # fig.canvas.draw() + fig.canvas.draw() - if self.save: # and self.data.ended: - self.contador += 1 + if self.save and self.data.ended: channels = range(self.nrows) if self.oneFigure: label = '' @@ -445,12 +445,11 @@ class PlotData(Operation, Process): label = '_{}'.format(channels[n]) figname = os.path.join( self.save, - '{}{}_{}{}.png'.format( + '{}{}_{}.png'.format( self.CODE, label, self.getDateTime(self.saveTime).strftime( - '%y%m%d_%H%M%S'), - str(self.contador), + '%y%m%d_%H%M%S'), ) ) log.log('Saving figure: {}'.format(figname), self.name) @@ -898,10 +897,11 @@ class PlotParamData(PlotRTIData): self.nplots += 1 self.ylabel = 'Height [Km]' - self.titles = self.data.parameters \ - if self.data.parameters else ['Param {}'.format(x) for x in xrange(self.nrows)] - if self.showSNR: - self.titles.append('SNR') + if not self.titles: + self.titles = self.data.parameters \ + if self.data.parameters else ['Param {}'.format(x) for x in xrange(self.nrows)] + if self.showSNR: + self.titles.append('SNR') def plot(self): self.data.normalize_heights() diff --git a/schainpy/model/io/__init__.py b/schainpy/model/io/__init__.py index dc9a6cd..f125f45 100644 --- a/schainpy/model/io/__init__.py +++ b/schainpy/model/io/__init__.py @@ -18,3 +18,4 @@ from jroIO_madrigal import * from bltrIO_param import * from jroIO_bltr import * from jroIO_mira35c import * +from julIO_param import * \ No newline at end of file diff --git a/schainpy/model/io/julIO_param.py b/schainpy/model/io/julIO_param.py new file mode 100644 index 0000000..208c0ba --- /dev/null +++ b/schainpy/model/io/julIO_param.py @@ -0,0 +1,343 @@ +''' +Created on Set 10, 2017 + +@author: Juan C. Espinoza +''' + + +import os +import sys +import time +import glob +import datetime + +import numpy + +from schainpy.model.proc.jroproc_base import ProcessingUnit +from schainpy.model.data.jrodata import Parameters +from schainpy.model.io.jroIO_base import JRODataReader, isNumber +from schainpy.utils import log + +FILE_HEADER_STRUCTURE = numpy.dtype([ + ('year', 'f'), + ('doy', 'f'), + ('nint', 'f'), + ('navg', 'f'), + ('fh', 'f'), + ('dh', 'f'), + ('nheights', 'f'), + ('ipp', 'f') +]) + +REC_HEADER_STRUCTURE = numpy.dtype([ + ('magic', 'f'), + ('hours', 'f'), + ('interval', 'f'), + ('h0', 'f'), + ('nheights', 'f'), + ('snr1', 'f'), + ('snr2', 'f'), + ('snr', 'f'), +]) + +DATA_STRUCTURE = numpy.dtype([ + ('range', ' dateFile) or (endDate < dateFile): + continue + + self.fileList.append(thisFile) + self.dateFileList.append(dateFile) + + return + + def setNextFile(self): + + file_id = self.fileIndex + + if file_id == len(self.fileList): + log.success('No more files in the folder', self.name) + self.flagNoMoreFiles = 1 + return 0 + + log.success('Opening {}'.format(self.fileList[file_id]), self.name) + filename = os.path.join(self.path, self.fileList[file_id]) + + dirname, name = os.path.split(filename) + self.siteFile = name.split('.')[0] + if self.filename is not None: + self.fp.close() + self.filename = filename + self.fp = open(self.filename, 'rb') + + self.header_file = numpy.fromfile(self.fp, FILE_HEADER_STRUCTURE, 1) + yy = self.header_file['year'] - 1900 * (self.header_file['year'] > 3000) + self.year = int(yy + 1900 * (yy < 1000)) + self.doy = int(self.header_file['doy']) + self.dH = round(self.header_file['dh'], 2) + self.ipp = round(self.header_file['ipp'], 2) + self.sizeOfFile = os.path.getsize(self.filename) + self.counter_records = 0 + self.flagIsNewFile = 0 + self.fileIndex += 1 + + return 1 + + def readNextBlock(self): + + while True: + if not self.readBlock(): + self.flagIsNewFile = 1 + if not self.setNextFile(): + return 0 + + if (self.datatime < datetime.datetime.combine(self.startDate, self.startTime)) or \ + (self.datatime > datetime.datetime.combine(self.endDate, self.endTime)): + log.warning( + 'Reading Record No. {} -> {} [Skipping]'.format( + self.counter_records, + self.datatime.ctime()), + self.name) + continue + break + + log.log('Reading Record No. {} -> {}'.format( + self.counter_records, + self.datatime.ctime()), self.name) + + return 1 + + def readBlock(self): + + pointer = self.fp.tell() + heights, dt = self.readHeader() + self.fp.seek(pointer) + buffer_h = [] + buffer_d = [] + while True: + pointer = self.fp.tell() + if pointer == self.sizeOfFile: + return 0 + heights, datatime = self.readHeader() + if dt == datatime: + buffer_h.append(heights) + buffer_d.append(self.readData(len(heights))) + continue + self.fp.seek(pointer) + break + + if dt.date() > self.datatime.date(): + self.flagDiscontinuousBlock = 1 + self.datatime = dt + self.time = (dt - datetime.datetime(1970, 1, 1)).total_seconds() + time.timezone + self.heights = numpy.concatenate(buffer_h) + self.buffer = numpy.zeros((5, len(self.heights))) + numpy.nan + self.buffer[0, :] = numpy.concatenate([buf[0] for buf in buffer_d]) + self.buffer[1, :] = numpy.concatenate([buf[1] for buf in buffer_d]) + self.buffer[2, :] = numpy.concatenate([buf[2] for buf in buffer_d]) + self.buffer[3, :] = numpy.concatenate([buf[3] for buf in buffer_d]) + self.buffer[4, :] = numpy.concatenate([buf[4] for buf in buffer_d]) + + self.counter_records += 1 + + return 1 + + def readHeader(self): + ''' + Parse recordHeader + ''' + + self.header_rec = numpy.fromfile(self.fp, REC_HEADER_STRUCTURE, 1) + self.interval = self.header_rec['interval'] + if self.header_rec['magic'] == 888.: + self.header_rec['h0'] = round(self.header_rec['h0'], 2) + nheights = int(self.header_rec['nheights']) + hours = float(self.header_rec['hours'][0]) + heights = numpy.arange(nheights) * self.dH + self.header_rec['h0'] + datatime = datetime.datetime(self.year, 1, 1) + datetime.timedelta(days=self.doy-1, hours=hours) + return heights, datatime + else: + return False + + def readData(self, N): + ''' + Parse data + ''' + + buffer = numpy.fromfile(self.fp, 'f', 8*N).reshape(N, 8) + + pow0 = buffer[:, 0] + pow1 = buffer[:, 1] + acf0 = (buffer[:,2] + buffer[:,3]*1j) / pow0 + acf1 = (buffer[:,4] + buffer[:,5]*1j) / pow1 + dccf = (buffer[:,6] + buffer[:,7]*1j) / (pow0*pow1) + + ### SNR + sno = (pow0 + pow1 - self.header_rec['snr']) / self.header_rec['snr'] + sno10 = numpy.log10(sno) + # dsno = 1.0 / numpy.sqrt(self.header_file['nint'] * self.header_file['navg']) * (1 + (1 / sno)) + + ### Vertical Drift + sp = numpy.sqrt(numpy.abs(acf0)*numpy.abs(acf1)) + sp[numpy.where(numpy.abs(sp) >= 1.0)] = numpy.sqrt(0.9999) + + vzo = -numpy.arctan2(acf0.imag + acf1.imag,acf0.real + acf1.real)*1.5E5*1.5/(self.ipp*numpy.pi) + dvzo = numpy.sqrt(1.0 - sp*sp)*0.338*1.5E5/(numpy.sqrt(self.header_file['nint']*self.header_file['navg'])*sp*self.ipp) + err = numpy.where(dvzo <= 0.1) + dvzo[err] = 0.1 + + #Zonal Drifts + dt = self.header_file['nint']*self.ipp / 1.5E5 + coh = numpy.sqrt(numpy.abs(dccf)) + err = numpy.where(coh >= 1.0) + coh[err] = numpy.sqrt(0.99999) + + err = numpy.where(coh <= 0.1) + coh[err] = numpy.sqrt(0.1) + + vxo = numpy.arctan2(dccf.imag, dccf.real)*self.header_rec['h0']*1.0E3/(self.kd*dt) + dvxo = numpy.sqrt(1.0 - coh*coh)*self.header_rec['h0']*1.0E3/(numpy.sqrt(self.header_file['nint']*self.header_file['navg'])*coh*self.kd*dt) + + err = numpy.where(dvxo <= 0.1) + dvxo[err] = 0.1 + + return vzo, dvzo, vxo, dvxo, sno10 + + def set_output(self): + ''' + Storing data from databuffer to dataOut object + ''' + + self.dataOut.data_SNR = self.buffer[4].reshape(1, -1) + self.dataOut.heightList = self.heights + self.dataOut.data_param = self.buffer[0:4,] + self.dataOut.utctimeInit = self.time + self.dataOut.utctime = self.time + self.dataOut.useLocalTime = True + self.dataOut.paramInterval = self.interval + self.dataOut.timezone = self.timezone + self.dataOut.sizeOfFile = self.sizeOfFile + self.dataOut.flagNoData = False + self.dataOut.flagDiscontinuousBlock = self.flagDiscontinuousBlock + + def getData(self): + ''' + Storing data from databuffer to dataOut object + ''' + if self.flagNoMoreFiles: + self.dataOut.flagNoData = True + log.success('No file left to process', self.name) + return 0 + + if not self.readNextBlock(): + self.dataOut.flagNoData = True + return 0 + + self.set_output() + + return 1 diff --git a/schainpy/model/utils/jroutils_publish.py b/schainpy/model/utils/jroutils_publish.py index 56d6163..a1236e6 100644 --- a/schainpy/model/utils/jroutils_publish.py +++ b/schainpy/model/utils/jroutils_publish.py @@ -138,12 +138,11 @@ class Data(object): return self.data[key][self.__times[0]].shape return (0,) - def update(self, dataOut): + def update(self, dataOut, tm): ''' Update data object with new dataOut ''' - - tm = dataOut.utctime + if tm in self.__times: return @@ -578,7 +577,7 @@ class PlotterReceiver(ProcessingUnit, Process): while True: dataOut = self.receiver.recv_pyobj() if not dataOut.flagNoData: - if dataOut.type == 'Parameters': + if dataOut.type == 'Parameters': tm = dataOut.utctimeInit else: tm = dataOut.utctime @@ -599,7 +598,7 @@ class PlotterReceiver(ProcessingUnit, Process): self.data.setup() self.dates.append(dt) - self.data.update(dataOut) + self.data.update(dataOut, tm) if dataOut.finished is True: self.connections -= 1