##// END OF EJS Templates
sin matplotlib de los modulos de claire
José Chávez -
r1075:eb67833bff79
parent child
Show More
@@ -5,8 +5,8
5 5 # from schainpy.model.utils.jroutils import *
6 6 # from schainpy.serializer import *
7 7
8 from graphics import *
8 9 from data import *
9 10 from io import *
10 11 from proc import *
11 from graphics import *
12 12 from utils import *
@@ -4,7 +4,7 import sys
4 4 import matplotlib
5 5
6 6 if 'linux' in sys.platform:
7 matplotlib.use("TKAgg")
7 matplotlib.use("TkAgg")
8 8
9 9 if 'darwin' in sys.platform:
10 10 matplotlib.use('TKAgg')
@@ -21,15 +21,19 from matplotlib.ticker import FuncFormatter, LinearLocator
21 21 # create jro colormap
22 22
23 23 jet_values = matplotlib.pyplot.get_cmap("jet", 100)(numpy.arange(100))[10:90]
24 blu_values = matplotlib.pyplot.get_cmap("seismic_r", 20)(numpy.arange(20))[10:15]
25 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list("jro", numpy.vstack((blu_values, jet_values)))
24 blu_values = matplotlib.pyplot.get_cmap(
25 "seismic_r", 20)(numpy.arange(20))[10:15]
26 ncmap = matplotlib.colors.LinearSegmentedColormap.from_list(
27 "jro", numpy.vstack((blu_values, jet_values)))
26 28 matplotlib.pyplot.register_cmap(cmap=ncmap)
27 29
30
28 31 def createFigure(id, wintitle, width, height, facecolor="w", show=True, dpi = 80):
29 32
30 33 matplotlib.pyplot.ioff()
31 34
32 fig = matplotlib.pyplot.figure(num=id, facecolor=facecolor, figsize=(1.0*width/dpi, 1.0*height/dpi))
35 fig = matplotlib.pyplot.figure(num=id, facecolor=facecolor, figsize=(
36 1.0 * width / dpi, 1.0 * height / dpi))
33 37 fig.canvas.manager.set_window_title(wintitle)
34 38 # fig.canvas.manager.resize(width, height)
35 39 matplotlib.pyplot.ion()
@@ -39,6 +43,7 def createFigure(id, wintitle, width, height, facecolor="w", show=True, dpi = 80
39 43
40 44 return fig
41 45
46
42 47 def closeFigure(show=False, fig=None):
43 48
44 49 # matplotlib.pyplot.ioff()
@@ -60,24 +65,29 def closeFigure(show=False, fig=None):
60 65
61 66 return
62 67
68
63 69 def saveFigure(fig, filename):
64 70
65 71 # matplotlib.pyplot.ioff()
66 72 fig.savefig(filename, dpi=matplotlib.pyplot.gcf().dpi)
67 73 # matplotlib.pyplot.ion()
68 74
75
69 76 def clearFigure(fig):
70 77
71 78 fig.clf()
72 79
80
73 81 def setWinTitle(fig, title):
74 82
75 83 fig.canvas.manager.set_window_title(title)
76 84
85
77 86 def setTitle(fig, title):
78 87
79 88 fig.suptitle(title)
80 89
90
81 91 def createAxes(fig, nrow, ncol, xpos, ypos, colspan, rowspan, polar=False):
82 92
83 93 matplotlib.pyplot.ioff()
@@ -91,6 +101,7 def createAxes(fig, nrow, ncol, xpos, ypos, colspan, rowspan, polar=False):
91 101 matplotlib.pyplot.ion()
92 102 return axes
93 103
104
94 105 def setAxesText(ax, text):
95 106
96 107 ax.annotate(text,
@@ -100,17 +111,18 def setAxesText(ax, text):
100 111 verticalalignment = 'top',
101 112 fontsize = 10)
102 113
114
103 115 def printLabels(ax, xlabel, ylabel, title):
104 116
105 117 ax.set_xlabel(xlabel, size=11)
106 118 ax.set_ylabel(ylabel, size=11)
107 119 ax.set_title(title, size=8)
108 120
121
109 122 def createPline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='',
110 123 ticksize=9, xtick_visible=True, ytick_visible=True,
111 124 nxticks=4, nyticks=10,
112 125 grid=None,color='blue'):
113
114 126 """
115 127
116 128 Input:
@@ -130,7 +142,8 def createPline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title=''
130 142 xtickspos = numpy.array([float("%.1f"%i) for i in xtickspos])
131 143 ax.set_xticks(xtickspos)
132 144 else:
133 xtickspos = numpy.arange(nxticks)*int((xmax-xmin)/(nxticks)) + int(xmin)
145 xtickspos = numpy.arange(nxticks) * \
146 int((xmax - xmin) / (nxticks)) + int(xmin)
134 147 # xtickspos = numpy.arange(nxticks)*float(xmax-xmin)/float(nxticks) + int(xmin)
135 148 ax.set_xticks(xtickspos)
136 149
@@ -168,10 +181,12 def createPline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title=''
168 181
169 182 return iplot
170 183
184
171 185 def set_linedata(ax, x, y, idline):
172 186
173 187 ax.lines[idline].set_data(x,y)
174 188
189
175 190 def pline(iplot, x, y, xlabel='', ylabel='', title=''):
176 191
177 192 ax = iplot.axes
@@ -180,6 +195,7 def pline(iplot, x, y, xlabel='', ylabel='', title=''):
180 195
181 196 set_linedata(ax, x, y, idline=0)
182 197
198
183 199 def addpline(ax, x, y, color, linestyle, lw):
184 200
185 201 ax.plot(x,y,color=color,linestyle=linestyle,lw=lw)
@@ -235,13 +251,15 def createPcolor(ax, x, y, z, xmin, xmax, ymin, ymax, zmin, zmax,
235 251
236 252 if XAxisAsTime:
237 253
238 func = lambda x, pos: ('%s') %(datetime.datetime.utcfromtimestamp(x).strftime("%H:%M:%S"))
254 def func(x, pos): return ('%s') % (
255 datetime.datetime.utcfromtimestamp(x).strftime("%H:%M:%S"))
239 256 ax.xaxis.set_major_formatter(FuncFormatter(func))
240 257 ax.xaxis.set_major_locator(LinearLocator(7))
241 258
242 259 matplotlib.pyplot.ion()
243 260 return imesh
244 261
262
245 263 def pcolor(imesh, z, xlabel='', ylabel='', title=''):
246 264
247 265 z = z.T
@@ -249,11 +267,14 def pcolor(imesh, z, xlabel='', ylabel='', title=''):
249 267 printLabels(ax, xlabel, ylabel, title)
250 268 imesh.set_array(z.ravel())
251 269
270
252 271 def addpcolor(ax, x, y, z, zmin, zmax, xlabel='', ylabel='', title='', colormap='jet'):
253 272
254 273 printLabels(ax, xlabel, ylabel, title)
255 274
256 ax.pcolormesh(x,y,z.T,vmin=zmin,vmax=zmax, cmap=matplotlib.pyplot.get_cmap(colormap))
275 ax.pcolormesh(x, y, z.T, vmin=zmin, vmax=zmax,
276 cmap=matplotlib.pyplot.get_cmap(colormap))
277
257 278
258 279 def addpcolorbuffer(ax, x, y, z, zmin, zmax, xlabel='', ylabel='', title='', colormap='jet'):
259 280
@@ -266,14 +287,13 def addpcolorbuffer(ax, x, y, z, zmin, zmax, xlabel='', ylabel='', title='', col
266 287 cmap=matplotlib.pyplot.get_cmap(colormap)
267 288 cmap.set_bad('black', 1.)
268 289
269
270 290 ax.pcolormesh(x,y,z.T,vmin=zmin,vmax=zmax, cmap=cmap)
271 291
292
272 293 def createPmultiline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
273 294 ticksize=9, xtick_visible=True, ytick_visible=True,
274 295 nxticks=4, nyticks=10,
275 296 grid=None):
276
277 297 """
278 298
279 299 Input:
@@ -289,7 +309,8 def createPmultiline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', tit
289 309 ax.set_ylim([ymin,ymax])
290 310 printLabels(ax, xlabel, ylabel, title)
291 311
292 xtickspos = numpy.arange(nxticks)*int((xmax-xmin)/(nxticks)) + int(xmin)
312 xtickspos = numpy.arange(nxticks) * \
313 int((xmax - xmin) / (nxticks)) + int(xmin)
293 314 ax.set_xticks(xtickspos)
294 315
295 316 for tick in ax.get_xticklabels():
@@ -334,11 +355,11 def pmultiline(iplot, x, y, xlabel='', ylabel='', title=''):
334 355 line = ax.lines[i]
335 356 line.set_data(x[i,:],y)
336 357
358
337 359 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
338 360 ticksize=9, xtick_visible=True, ytick_visible=True,
339 361 nxticks=4, nyticks=10, marker='.', markersize=10, linestyle="None",
340 362 grid=None, XAxisAsTime=False):
341
342 363 """
343 364
344 365 Input:
@@ -355,7 +376,8 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel=''
355 376 leg = ax.legend(lines, legendlabels,
356 377 loc='upper right', bbox_to_anchor=(1.16, 1), borderaxespad=0)
357 378
358 for label in leg.get_texts(): label.set_fontsize(9)
379 for label in leg.get_texts():
380 label.set_fontsize(9)
359 381
360 382 ax.set_xlim([xmin,xmax])
361 383 ax.set_ylim([ymin,ymax])
@@ -393,7 +415,8 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel=''
393 415
394 416 if XAxisAsTime:
395 417
396 func = lambda x, pos: ('%s') %(datetime.datetime.utcfromtimestamp(x).strftime("%H:%M:%S"))
418 def func(x, pos): return ('%s') % (
419 datetime.datetime.utcfromtimestamp(x).strftime("%H:%M:%S"))
397 420 ax.xaxis.set_major_formatter(FuncFormatter(func))
398 421 ax.xaxis.set_major_locator(LinearLocator(7))
399 422
@@ -401,6 +424,7 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel=''
401 424
402 425 return iplot
403 426
427
404 428 def pmultilineyaxis(iplot, x, y, xlabel='', ylabel='', title=''):
405 429
406 430 ax = iplot.axes
@@ -411,6 +435,7 def pmultilineyaxis(iplot, x, y, xlabel='', ylabel='', title=''):
411 435 line = ax.lines[i]
412 436 line.set_data(x,y[i,:])
413 437
438
414 439 def createPolar(ax, x, y,
415 440 xlabel='', ylabel='', title='', ticksize = 9,
416 441 colormap='jet',cblabel='', cbsize="5%",
@@ -444,9 +469,9 def createPolar(ax, x, y,
444 469
445 470 matplotlib.pyplot.ion()
446 471
447
448 472 return iplot
449 473
474
450 475 def polar(iplot, x, y, xlabel='', ylabel='', title=''):
451 476
452 477 ax = iplot.axes
@@ -456,6 +481,7 def polar(iplot, x, y, xlabel='', ylabel='', title=''):
456 481
457 482 set_linedata(ax, x, y, idline=0)
458 483
484
459 485 def draw(fig):
460 486
461 487 if type(fig) == 'int':
@@ -463,6 +489,7 def draw(fig):
463 489
464 490 fig.canvas.draw()
465 491
492
466 493 def pause(interval=0.000001):
467 494
468 495 matplotlib.pyplot.pause(interval)
@@ -1,4 +1,5
1 import os, sys
1 import os
2 import sys
2 3 import glob
3 4 import fnmatch
4 5 import datetime
@@ -6,9 +7,7 import time
6 7 import re
7 8 import h5py
8 9 import numpy
9 import matplotlib.pyplot as plt
10 10
11 import pylab as plb
12 11 from scipy.optimize import curve_fit
13 12 from scipy import asarray as ar,exp
14 13 from scipy import stats
@@ -31,15 +30,20 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
31 30 from numpy import imag, shape, NaN
32 31
33 32
34 startFp = open('/home/erick/Documents/MIRA35C/20160117/20160117_0000.zspc',"rb")
33 startFp = open(
34 '/home/erick/Documents/MIRA35C/20160117/20160117_0000.zspc', "rb")
35 35
36 36
37 37 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
38 38 ('Hname',numpy.str_,32), #Original file name
39 ('Htime',numpy.str_,32), #Date and time when the file was created
40 ('Hoper',numpy.str_,64), #Name of operator who created the file
41 ('Hplace',numpy.str_,128), #Place where the measurements was carried out
42 ('Hdescr',numpy.str_,256), #Description of measurements
39 # Date and time when the file was created
40 ('Htime', numpy.str_, 32),
41 # Name of operator who created the file
42 ('Hoper', numpy.str_, 64),
43 # Place where the measurements was carried out
44 ('Hplace', numpy.str_, 128),
45 # Description of measurements
46 ('Hdescr', numpy.str_, 256),
43 47 ('Hdummy',numpy.str_,512), #Reserved space
44 48 #Main chunk
45 49 ('Msign','<i4'), #Main chunk signature FZKF or NUIG
@@ -50,19 +54,28 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
50 54 ('PPARprf','<i4'), #Pulse repetition frequency
51 55 ('PPARpdr','<i4'), #Pulse duration
52 56 ('PPARsft','<i4'), #FFT length
53 ('PPARavc','<i4'), #Number of spectral (in-coherent) averages
54 ('PPARihp','<i4'), #Number of lowest range gate for moment estimation
55 ('PPARchg','<i4'), #Count for gates for moment estimation
56 ('PPARpol','<i4'), #switch on/off polarimetric measurements. Should be 1.
57 # Number of spectral (in-coherent) averages
58 ('PPARavc', '<i4'),
59 # Number of lowest range gate for moment estimation
60 ('PPARihp', '<i4'),
61 # Count for gates for moment estimation
62 ('PPARchg', '<i4'),
63 # switch on/off polarimetric measurements. Should be 1.
64 ('PPARpol', '<i4'),
57 65 #Service DSP parameters
58 ('SPARatt','<i4'), #STC attenuation on the lowest ranges on/off
66 # STC attenuation on the lowest ranges on/off
67 ('SPARatt', '<i4'),
59 68 ('SPARtx','<i4'), #OBSOLETE
60 69 ('SPARaddGain0','<f4'), #OBSOLETE
61 70 ('SPARaddGain1','<f4'), #OBSOLETE
62 ('SPARwnd','<i4'), #Debug only. It normal mode it is 0.
63 ('SPARpos','<i4'), #Delay between sync pulse and tx pulse for phase corr, ns
64 ('SPARadd','<i4'), #"add to pulse" to compensate for delay between the leading edge of driver pulse and envelope of the RF signal.
65 ('SPARlen','<i4'), #Time for measuring txn pulse phase. OBSOLETE
71 # Debug only. It normal mode it is 0.
72 ('SPARwnd', '<i4'),
73 # Delay between sync pulse and tx pulse for phase corr, ns
74 ('SPARpos', '<i4'),
75 # "add to pulse" to compensate for delay between the leading edge of driver pulse and envelope of the RF signal.
76 ('SPARadd', '<i4'),
77 # Time for measuring txn pulse phase. OBSOLETE
78 ('SPARlen', '<i4'),
66 79 ('SPARcal','<i4'), #OBSOLETE
67 80 ('SPARnos','<i4'), #OBSOLETE
68 81 ('SPARof0','<i4'), #detection threshold
@@ -73,19 +86,23 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
73 86 ('SPARtst','<i4'), #OBSOLETE
74 87 ('SPARcor','<i4'), #OBSOLETE
75 88 ('SPARofs','<i4'), #OBSOLETE
76 ('SPARhsn','<i4'), #Hildebrand div noise detection on noise gate
77 ('SPARhsa','<f4'), #Hildebrand div noise detection on all gates
89 # Hildebrand div noise detection on noise gate
90 ('SPARhsn', '<i4'),
91 # Hildebrand div noise detection on all gates
92 ('SPARhsa', '<f4'),
78 93 ('SPARcalibPow_M','<f4'), #OBSOLETE
79 94 ('SPARcalibSNR_M','<f4'), #OBSOLETE
80 95 ('SPARcalibPow_S','<f4'), #OBSOLETE
81 96 ('SPARcalibSNR_S','<f4'), #OBSOLETE
82 ('SPARrawGate1','<i4'), #Lowest range gate for spectra saving Raw_Gate1 >=5
83 ('SPARrawGate2','<i4'), #Number of range gates with atmospheric signal
84 ('SPARraw','<i4'), #flag - IQ or spectra saving on/off
97 # Lowest range gate for spectra saving Raw_Gate1 >=5
98 ('SPARrawGate1', '<i4'),
99 # Number of range gates with atmospheric signal
100 ('SPARrawGate2', '<i4'),
101 # flag - IQ or spectra saving on/off
102 ('SPARraw', '<i4'),
85 103 ('SPARprc','<i4'),]) #flag - Moment estimation switched on/off
86 104
87 105
88
89 106 self.Hname= None
90 107 self.Htime= None
91 108 self.Hoper= None
@@ -136,7 +153,6 self.SPARraw=None
136 153 self.SPARprc=None
137 154
138 155
139
140 156 header = numpy.fromfile(fp, FILE_HEADER,1)
141 157 ''' numpy.fromfile(file, dtype, count, sep='')
142 158 file : file or str
@@ -207,9 +223,8 SPARraw=header['SPARraw'][0]
207 223 SPARprc=header['SPARprc'][0]
208 224
209 225
210
211 226 SRVI_STRUCTURE = numpy.dtype([
212 ('frame_cnt','<u4'),#
227 ('frame_cnt', '<u4'),
213 228 ('time_t','<u4'), #
214 229 ('tpow','<f4'), #
215 230 ('npw1','<f4'), #
@@ -225,19 +240,18 SRVI_STRUCTURE = numpy.dtype([
225 240 ('azivel','<f4'), #
226 241 ('elvpos','<f4'), #
227 242 ('elvvel','<f4'), #
228 ('northAngle','<f4'), #
243 ('northAngle', '<f4'),
229 244 ('microsec','<u4'), #
230 245 ('azisetvel','<f4'), #
231 246 ('elvsetpos','<f4'), #
232 247 ('RadarConst','<f4'),]) #
233 248
234 249 JUMP_STRUCTURE = numpy.dtype([
235 ('jump','<u140'),#
236 ('SizeOfDataBlock1',numpy.str_,32),#
237 ('jump','<i4'),#
238 ('DataBlockTitleSRVI1',numpy.str_,32),#
239 ('SizeOfSRVI1','<i4'),])#
240
250 ('jump', '<u140'),
251 ('SizeOfDataBlock1', numpy.str_, 32),
252 ('jump', '<i4'),
253 ('DataBlockTitleSRVI1', numpy.str_, 32),
254 ('SizeOfSRVI1', '<i4'), ])
241 255
242 256
243 257 #frame_cnt=0, time_t= 0, tpow=0, npw1=0, npw2=0,
@@ -269,7 +283,6 elvsetpos = elvsetpos
269 283 RadarConst5 = RadarConst
270 284
271 285
272
273 286 #print fp
274 287 #startFp = open('/home/erick/Documents/Data/huancayo.20161019.22.fdt',"rb") #The method tell() returns the current position of the file read/write pointer within the file.
275 288 #startFp = open(fp,"rb") #The method tell() returns the current position of the file read/write pointer within the file.
@@ -283,7 +296,7 print 'Posicion del bloque: ',OffRHeader
283 296
284 297 header = numpy.fromfile(startFp,SRVI_STRUCTURE,1)
285 298
286 self.frame_cnt = header['frame_cnt'][0]#
299 self.frame_cnt = header['frame_cnt'][0]
287 300 self.time_t = header['frame_cnt'][0] #
288 301 self.tpow = header['frame_cnt'][0] #
289 302 self.npw1 = header['frame_cnt'][0] #
@@ -316,6 +329,3 endFp = self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
316 329 print '=============================================='
317 330
318 331 print '=============================================='
319
320
321 No newline at end of file
@@ -1,4 +1,5
1 import os, sys
1 import os
2 import sys
2 3 import glob
3 4 import fnmatch
4 5 import datetime
@@ -6,7 +7,6 import time
6 7 import re
7 8 import h5py
8 9 import numpy
9 import matplotlib.pyplot as plt
10 10
11 11 import pylab as plb
12 12 from scipy.optimize import curve_fit
@@ -37,7 +37,6 class Header(object):
37 37 def __init__(self):
38 38 raise NotImplementedError
39 39
40
41 40 def read(self):
42 41
43 42 raise NotImplementedError
@@ -67,23 +66,23 class Header(object):
67 66 #print message
68 67
69 68
70
71
72
73 69 FILE_STRUCTURE = numpy.dtype([ #HEADER 48bytes
74 70 ('FileMgcNumber','<u4'), #0x23020100
75 ('nFDTdataRecors','<u4'), #No Of FDT data records in this file (0 or more)
71 # No Of FDT data records in this file (0 or more)
72 ('nFDTdataRecors', '<u4'),
76 73 ('OffsetStartHeader','<u4'),
77 74 ('RadarUnitId','<u4'),
78 75 ('SiteName',numpy.str_,32), #Null terminated
79 76 ])
80 77
78
81 79 class FileHeaderBLTR(Header):
82 80
83 81 def __init__(self):
84 82
85 83 self.FileMgcNumber= 0 #0x23020100
86 self.nFDTdataRecors=0 #No Of FDT data records in this file (0 or more)
84 # No Of FDT data records in this file (0 or more)
85 self.nFDTdataRecors = 0
87 86 self.RadarUnitId= 0
88 87 self.OffsetStartHeader=0
89 88 self.SiteName= ""
@@ -99,7 +98,6 class FileHeaderBLTR(Header):
99 98 print 'puntero file header', startFp.tell()
100 99 print ' '
101 100
102
103 101 ''' numpy.fromfile(file, dtype, count, sep='')
104 102 file : file or str
105 103 Open file object or filename.
@@ -119,23 +117,20 class FileHeaderBLTR(Header):
119 117
120 118 '''
121 119
122
123
124 120 self.FileMgcNumber= hex(header['FileMgcNumber'][0])
125 self.nFDTdataRecors=int(header['nFDTdataRecors'][0]) #No Of FDT data records in this file (0 or more)
121 # No Of FDT data records in this file (0 or more)
122 self.nFDTdataRecors = int(header['nFDTdataRecors'][0])
126 123 self.RadarUnitId= int(header['RadarUnitId'][0])
127 124 self.OffsetStartHeader= int(header['OffsetStartHeader'][0])
128 125 self.SiteName= str(header['SiteName'][0])
129 126
130 127 #print 'Numero de bloques', self.nFDTdataRecors
131 128
132
133 129 if self.size <48:
134 130 return 0
135 131
136 132 return 1
137 133
138
139 134 def write(self, fp):
140 135
141 136 headerTuple = (self.FileMgcNumber,
@@ -144,7 +139,6 class FileHeaderBLTR(Header):
144 139 self.SiteName,
145 140 self.size)
146 141
147
148 142 header = numpy.array(headerTuple, FILE_STRUCTURE)
149 143 # numpy.array(object, dtype=None, copy=True, order=None, subok=False, ndmin=0)
150 144 header.tofile(fp)
@@ -166,56 +160,92 class FileHeaderBLTR(Header):
166 160 return 1
167 161
168 162
169
170
171
172 163 RECORD_STRUCTURE = numpy.dtype([ #RECORD HEADER 180+20N bytes
173 164 ('RecMgcNumber','<u4'), #0x23030001
174 165 ('RecCounter','<u4'), #Record counter(0,1, ...)
175 ('Off2StartNxtRec','<u4'), #Offset to start of next record form start of this record
176 ('Off2StartData','<u4'), #Offset to start of data from start of this record
177 ('nUtime','<i4'), #Epoch time stamp of start of acquisition (seconds)
178 ('nMilisec','<u4'), #Millisecond component of time stamp (0,...,999)
179 ('ExpTagName',numpy.str_,32), #Experiment tag name (null terminated)
180 ('ExpComment',numpy.str_,32), #Experiment comment (null terminated)
181 ('SiteLatDegrees','<f4'), #Site latitude (from GPS) in degrees (positive implies North)
182 ('SiteLongDegrees','<f4'), #Site longitude (from GPS) in degrees (positive implies East)
183 ('RTCgpsStatus','<u4'), #RTC GPS engine status (0=SEEK, 1=LOCK, 2=NOT FITTED, 3=UNAVAILABLE)
166 # Offset to start of next record form start of this record
167 ('Off2StartNxtRec', '<u4'),
168 # Offset to start of data from start of this record
169 ('Off2StartData', '<u4'),
170 # Epoch time stamp of start of acquisition (seconds)
171 ('nUtime', '<i4'),
172 # Millisecond component of time stamp (0,...,999)
173 ('nMilisec', '<u4'),
174 # Experiment tag name (null terminated)
175 ('ExpTagName', numpy.str_, 32),
176 # Experiment comment (null terminated)
177 ('ExpComment', numpy.str_, 32),
178 # Site latitude (from GPS) in degrees (positive implies North)
179 ('SiteLatDegrees', '<f4'),
180 # Site longitude (from GPS) in degrees (positive implies East)
181 ('SiteLongDegrees', '<f4'),
182 # RTC GPS engine status (0=SEEK, 1=LOCK, 2=NOT FITTED, 3=UNAVAILABLE)
183 ('RTCgpsStatus', '<u4'),
184 184 ('TransmitFrec','<u4'), #Transmit frequency (Hz)
185 185 ('ReceiveFrec','<u4'), #Receive frequency
186 ('FirstOsciFrec','<u4'), #First local oscillator frequency (Hz)
187 ('Polarisation','<u4'), #(0="O", 1="E", 2="linear 1", 3="linear2")
188 ('ReceiverFiltSett','<u4'), #Receiver filter settings (0,1,2,3)
189 ('nModesInUse','<u4'), #Number of modes in use (1 or 2)
190 ('DualModeIndex','<u4'), #Dual Mode index number for these data (0 or 1)
191 ('DualModeRange','<u4'), #Dual Mode range correction for these data (m)
192 ('nDigChannels','<u4'), #Number of digital channels acquired (2*N)
193 ('SampResolution','<u4'), #Sampling resolution (meters)
194 ('nHeights','<u4'), #Number of range gates sampled
195 ('StartRangeSamp','<u4'), #Start range of sampling (meters)
186 # First local oscillator frequency (Hz)
187 ('FirstOsciFrec', '<u4'),
188 # (0="O", 1="E", 2="linear 1", 3="linear2")
189 ('Polarisation', '<u4'),
190 # Receiver filter settings (0,1,2,3)
191 ('ReceiverFiltSett', '<u4'),
192 # Number of modes in use (1 or 2)
193 ('nModesInUse', '<u4'),
194 # Dual Mode index number for these data (0 or 1)
195 ('DualModeIndex', '<u4'),
196 # Dual Mode range correction for these data (m)
197 ('DualModeRange', '<u4'),
198 # Number of digital channels acquired (2*N)
199 ('nDigChannels', '<u4'),
200 # Sampling resolution (meters)
201 ('SampResolution', '<u4'),
202 # Number of range gates sampled
203 ('nHeights', '<u4'),
204 # Start range of sampling (meters)
205 ('StartRangeSamp', '<u4'),
196 206 ('PRFhz','<u4'), #PRF (Hz)
197 207 ('nCohInt','<u4'), #Integrations
198 ('nProfiles','<u4'), #Number of data points transformed
199 ('nChannels','<u4'), #Number of receive beams stored in file (1 or N)
208 # Number of data points transformed
209 ('nProfiles', '<u4'),
210 # Number of receive beams stored in file (1 or N)
211 ('nChannels', '<u4'),
200 212 ('nIncohInt','<u4'), #Number of spectral averages
201 ('FFTwindowingInd','<u4'), #FFT windowing index (0 = no window)
202 ('BeamAngleAzim','<f4'), #Beam steer angle (azimuth) in degrees (clockwise from true North)
203 ('BeamAngleZen','<f4'), #Beam steer angle (zenith) in degrees (0=> vertical)
204 ('AntennaCoord0','<f4'), #Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
205 ('AntennaAngl0','<f4'), #Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
206 ('AntennaCoord1','<f4'), #Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
207 ('AntennaAngl1','<f4'), #Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
208 ('AntennaCoord2','<f4'), #Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
209 ('AntennaAngl2','<f4'), #Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
210 ('RecPhaseCalibr0','<f4'), #Receiver phase calibration (degrees) - N values
211 ('RecPhaseCalibr1','<f4'), #Receiver phase calibration (degrees) - N values
212 ('RecPhaseCalibr2','<f4'), #Receiver phase calibration (degrees) - N values
213 ('RecAmpCalibr0','<f4'), #Receiver amplitude calibration (ratio relative to receiver one) - N values
214 ('RecAmpCalibr1','<f4'), #Receiver amplitude calibration (ratio relative to receiver one) - N values
215 ('RecAmpCalibr2','<f4'), #Receiver amplitude calibration (ratio relative to receiver one) - N values
216 ('ReceiverGaindB0','<i4'), #Receiver gains in dB - N values
217 ('ReceiverGaindB1','<i4'), #Receiver gains in dB - N values
218 ('ReceiverGaindB2','<i4'), #Receiver gains in dB - N values
213 # FFT windowing index (0 = no window)
214 ('FFTwindowingInd', '<u4'),
215 # Beam steer angle (azimuth) in degrees (clockwise from true North)
216 ('BeamAngleAzim', '<f4'),
217 # Beam steer angle (zenith) in degrees (0=> vertical)
218 ('BeamAngleZen', '<f4'),
219 # Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
220 ('AntennaCoord0', '<f4'),
221 # Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
222 ('AntennaAngl0', '<f4'),
223 # Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
224 ('AntennaCoord1', '<f4'),
225 # Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
226 ('AntennaAngl1', '<f4'),
227 # Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
228 ('AntennaCoord2', '<f4'),
229 # Antenna coordinates (Range(meters), Bearing(degrees)) - N pairs
230 ('AntennaAngl2', '<f4'),
231 # Receiver phase calibration (degrees) - N values
232 ('RecPhaseCalibr0', '<f4'),
233 # Receiver phase calibration (degrees) - N values
234 ('RecPhaseCalibr1', '<f4'),
235 # Receiver phase calibration (degrees) - N values
236 ('RecPhaseCalibr2', '<f4'),
237 # Receiver amplitude calibration (ratio relative to receiver one) - N values
238 ('RecAmpCalibr0', '<f4'),
239 # Receiver amplitude calibration (ratio relative to receiver one) - N values
240 ('RecAmpCalibr1', '<f4'),
241 # Receiver amplitude calibration (ratio relative to receiver one) - N values
242 ('RecAmpCalibr2', '<f4'),
243 # Receiver gains in dB - N values
244 ('ReceiverGaindB0', '<i4'),
245 # Receiver gains in dB - N values
246 ('ReceiverGaindB1', '<i4'),
247 # Receiver gains in dB - N values
248 ('ReceiverGaindB2', '<i4'),
219 249 ])
220 250
221 251
@@ -285,12 +315,11 class RecordHeaderBLTR(Header):
285 315 self.ReceiverGaindB2 = ReceiverGaindB2
286 316 self.OffsetStartHeader = 48
287 317
288
289
290 318 def RHread(self, fp):
291 319 #print fp
292 320 #startFp = open('/home/erick/Documents/Data/huancayo.20161019.22.fdt',"rb") #The method tell() returns the current position of the file read/write pointer within the file.
293 startFp = open(fp,"rb") #The method tell() returns the current position of the file read/write pointer within the file.
321 # The method tell() returns the current position of the file read/write pointer within the file.
322 startFp = open(fp, "rb")
294 323 #RecCounter=0
295 324 #Off2StartNxtRec=811248
296 325 OffRHeader= self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
@@ -298,7 +327,6 class RecordHeaderBLTR(Header):
298 327 print 'puntero Record Header', startFp.tell()
299 328 print ' '
300 329
301
302 330 startFp.seek(OffRHeader, os.SEEK_SET)
303 331
304 332 print ' '
@@ -416,11 +444,13 class RecordHeaderBLTR(Header):
416 444 print '=============================================='
417 445
418 446 if OffRHeader > endFp:
419 sys.stderr.write("Warning %s: Size value read from System Header is lower than it has to be\n" %fp)
447 sys.stderr.write(
448 "Warning %s: Size value read from System Header is lower than it has to be\n" % fp)
420 449 return 0
421 450
422 451 if OffRHeader < endFp:
423 sys.stderr.write("Warning %s: Size value read from System Header size is greater than it has to be\n" %fp)
452 sys.stderr.write(
453 "Warning %s: Size value read from System Header size is greater than it has to be\n" % fp)
424 454 return 0
425 455
426 456 return 1
@@ -436,7 +466,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
436 466 walk = None
437 467 isConfig = False
438 468
439
440 469 fileList= None
441 470
442 471 #metadata
@@ -448,8 +477,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
448 477 data= None
449 478 utctime= None
450 479
451
452
453 480 def __init__(self, **kwargs):
454 481
455 482 #Eliminar de la base la herencia
@@ -494,8 +521,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
494 521 self.dataOut.velocityY=[]
495 522 self.dataOut.velocityV=[]
496 523
497
498
499 524 def Files2Read(self, fp):
500 525 '''
501 526 Function that indicates the number of .fdt files that exist in the folder to be read.
@@ -503,8 +528,10 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
503 528 '''
504 529 #self.__checkPath()
505 530
506 ListaData=os.listdir(fp) #Gets the list of files within the fp address
507 ListaData=sorted(ListaData) #Sort the list of files from least to largest by names
531 # Gets the list of files within the fp address
532 ListaData = os.listdir(fp)
533 # Sort the list of files from least to largest by names
534 ListaData = sorted(ListaData)
508 535 nFiles=0 #File Counter
509 536 FileList=[] #A list is created that will contain the .fdt files
510 537 for IndexFile in ListaData :
@@ -517,7 +544,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
517 544
518 545 self.filenameList=FileList #List of files from least to largest by names
519 546
520
521 547 def run(self, **kwargs):
522 548 '''
523 549 This method will be the one that will initiate the data entry, will be called constantly.
@@ -531,7 +557,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
531 557 self.getData()
532 558 #print 'running'
533 559
534
535 560 def setup(self, path=None,
536 561 startDate=None,
537 562 endDate=None,
@@ -556,7 +581,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
556 581
557 582 pass
558 583
559
560 584 def getData(self):
561 585 '''
562 586 Before starting this function, you should check that there is still an unread file,
@@ -583,7 +607,6 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
583 607 #self.removeDC()
584 608 return self.dataOut.data_spc
585 609
586
587 610 def readFile(self,fp):
588 611 '''
589 612 You must indicate if you are reading in Online or Offline mode and load the
@@ -600,7 +623,8 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
600 623
601 624 if self.fileSelector < len(self.filenameList):
602 625
603 self.fpFile=str(fp)+'/'+str(self.filenameList[self.fileSelector])
626 self.fpFile = str(fp) + '/' + \
627 str(self.filenameList[self.fileSelector])
604 628 #print self.fpFile
605 629 fheader = FileHeaderBLTR()
606 630 fheader.FHread(self.fpFile) #Bltr FileHeader Reading
@@ -615,12 +639,15 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
615 639
616 640 def getVelRange(self, extrapoints=0):
617 641 Lambda= SPEED_OF_LIGHT/50000000
618 PRF = self.dataOut.PRF#1./(self.dataOut.ippSeconds * self.dataOut.nCohInt)
642 # 1./(self.dataOut.ippSeconds * self.dataOut.nCohInt)
643 PRF = self.dataOut.PRF
619 644 Vmax=-Lambda/(4.*(1./PRF)*self.dataOut.nCohInt*2.)
620 645 deltafreq = PRF / (self.nProfiles)
621 646 deltavel = (Vmax*2) / (self.nProfiles)
622 freqrange = deltafreq*(numpy.arange(self.nProfiles)-self.nProfiles/2.) - deltafreq/2
623 velrange = deltavel*(numpy.arange(self.nProfiles)-self.nProfiles/2.)
647 freqrange = deltafreq * \
648 (numpy.arange(self.nProfiles) - self.nProfiles / 2.) - deltafreq / 2
649 velrange = deltavel * \
650 (numpy.arange(self.nProfiles) - self.nProfiles / 2.)
624 651 return velrange
625 652
626 653 def readBlock(self):
@@ -661,7 +688,8 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
661 688
662 689 self.__firstHeigth=rheader.StartRangeSamp
663 690 self.__deltaHeigth=rheader.SampResolution
664 self.dataOut.heightList= self.__firstHeigth + numpy.array(range(self.nHeights))*self.__deltaHeigth
691 self.dataOut.heightList = self.__firstHeigth + \
692 numpy.array(range(self.nHeights)) * self.__deltaHeigth
665 693 self.dataOut.channelList = range(self.nChannels)
666 694 self.dataOut.nProfiles=rheader.nProfiles
667 695 self.dataOut.nIncohInt=rheader.nIncohInt
@@ -671,8 +699,10 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
671 699 self.dataOut.nFFTPoints=rheader.nProfiles
672 700 self.dataOut.utctime=rheader.nUtime
673 701 self.dataOut.timeZone=0
674 self.dataOut.normFactor= self.dataOut.nProfiles*self.dataOut.nIncohInt*self.dataOut.nCohInt
675 self.dataOut.outputInterval= self.dataOut.ippSeconds * self.dataOut.nCohInt * self.dataOut.nIncohInt * self.nProfiles
702 self.dataOut.normFactor = self.dataOut.nProfiles * \
703 self.dataOut.nIncohInt * self.dataOut.nCohInt
704 self.dataOut.outputInterval = self.dataOut.ippSeconds * \
705 self.dataOut.nCohInt * self.dataOut.nIncohInt * self.nProfiles
676 706
677 707 self.data_output=numpy.ones([3,rheader.nHeights])*numpy.NaN
678 708 print 'self.data_output', shape(self.data_output)
@@ -686,7 +716,8 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
686 716
687 717 #Procedure to take the pointer to where the date block starts
688 718 startDATA = open(self.fpFile,"rb")
689 OffDATA= self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec+self.Off2StartData
719 OffDATA = self.OffsetStartHeader + self.RecCounter * \
720 self.Off2StartNxtRec + self.Off2StartData
690 721 startDATA.seek(OffDATA, os.SEEK_SET)
691 722
692 723 def moving_average(x, N=2):
@@ -705,27 +736,26 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
705 736 y = rho * numpy.sin(phi)
706 737 return(x, y)
707 738
708
709
710
711 739 if self.DualModeIndex==self.ReadMode:
712 740
713 self.data_fft = numpy.fromfile( startDATA, [('complex','<c8')],self.nProfiles*self.nChannels*self.nHeights )
741 self.data_fft = numpy.fromfile(
742 startDATA, [('complex', '<c8')], self.nProfiles * self.nChannels * self.nHeights)
714 743
715 744 self.data_fft=self.data_fft.astype(numpy.dtype('complex'))
716 745
717 self.data_block=numpy.reshape(self.data_fft,(self.nHeights, self.nChannels, self.nProfiles ))
746 self.data_block = numpy.reshape(
747 self.data_fft, (self.nHeights, self.nChannels, self.nProfiles))
718 748
719 749 self.data_block = numpy.transpose(self.data_block, (1,2,0))
720 750
721 751 copy = self.data_block.copy()
722 752 spc = copy * numpy.conjugate(copy)
723 753
724 self.data_spc = numpy.absolute(spc) # valor absoluto o magnitud
754 self.data_spc = numpy.absolute(
755 spc) # valor absoluto o magnitud
725 756
726 757 factor = self.dataOut.normFactor
727 758
728
729 759 z = self.data_spc.copy()#/factor
730 760 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
731 761 #zdB = 10*numpy.log10(z)
@@ -737,13 +767,14 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
737 767
738 768 self.dataOut.data_spc=self.data_spc
739 769
740 self.noise = self.dataOut.getNoise(ymin_index=80, ymax_index=132)#/factor
770 self.noise = self.dataOut.getNoise(
771 ymin_index=80, ymax_index=132) # /factor
741 772 #noisedB = 10*numpy.log10(self.noise)
742 773
743
744 774 ySamples=numpy.ones([3,self.nProfiles])
745 775 phase=numpy.ones([3,self.nProfiles])
746 CSPCSamples=numpy.ones([3,self.nProfiles],dtype=numpy.complex_)
776 CSPCSamples = numpy.ones(
777 [3, self.nProfiles], dtype=numpy.complex_)
747 778 coherence=numpy.ones([3,self.nProfiles])
748 779 PhaseSlope=numpy.ones(3)
749 780 PhaseInter=numpy.ones(3)
@@ -766,13 +797,16 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
766 797 chan_index0 = self.dataOut.pairsList[i][0]
767 798 chan_index1 = self.dataOut.pairsList[i][1]
768 799
769 self.data_cspc[i,:,:]=cspc[chan_index0,:,:] * numpy.conjugate(cspc[chan_index1,:,:])
770
800 self.data_cspc[i, :, :] = cspc[chan_index0, :,
801 :] * numpy.conjugate(cspc[chan_index1, :, :])
771 802
772 803 '''Getting Eij and Nij'''
773 (AntennaX0,AntennaY0)=pol2cart(rheader.AntennaCoord0, rheader.AntennaAngl0*numpy.pi/180)
774 (AntennaX1,AntennaY1)=pol2cart(rheader.AntennaCoord1, rheader.AntennaAngl1*numpy.pi/180)
775 (AntennaX2,AntennaY2)=pol2cart(rheader.AntennaCoord2, rheader.AntennaAngl2*numpy.pi/180)
804 (AntennaX0, AntennaY0) = pol2cart(
805 rheader.AntennaCoord0, rheader.AntennaAngl0 * numpy.pi / 180)
806 (AntennaX1, AntennaY1) = pol2cart(
807 rheader.AntennaCoord1, rheader.AntennaAngl1 * numpy.pi / 180)
808 (AntennaX2, AntennaY2) = pol2cart(
809 rheader.AntennaCoord2, rheader.AntennaAngl2 * numpy.pi / 180)
776 810
777 811 E01=AntennaX0-AntennaX1
778 812 N01=AntennaY0-AntennaY1
@@ -783,7 +817,8 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
783 817 E12=AntennaX1-AntennaX2
784 818 N12=AntennaY1-AntennaY2
785 819
786 self.ChanDist= numpy.array([[E01, N01],[E02,N02],[E12,N12]])
820 self.ChanDist = numpy.array(
821 [[E01, N01], [E02, N02], [E12, N12]])
787 822
788 823 self.dataOut.ChanDist = self.ChanDist
789 824
@@ -1139,16 +1174,9 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1139 1174 #
1140 1175 # print ' '
1141 1176
1142
1143
1144 1177 self.BlockCounter+=2
1145 1178
1146 1179 else:
1147 1180 self.fileSelector+=1
1148 1181 self.BlockCounter=0
1149 1182 print "Next File"
1150
1151
1152
1153
1154
@@ -1,4 +1,5
1 import os, sys
1 import os
2 import sys
2 3 import glob
3 4 import fnmatch
4 5 import datetime
@@ -6,9 +7,7 import time
6 7 import re
7 8 import h5py
8 9 import numpy
9 import matplotlib.pyplot as plt
10 10
11 import pylab as plb
12 11 from scipy.optimize import curve_fit
13 12 from scipy import asarray as ar,exp
14 13 from scipy import stats
@@ -30,13 +29,11 from schainpy.model.proc.jroproc_base import ProcessingUnit, Operation
30 29 from numpy import imag, shape, NaN, empty
31 30
32 31
33
34 32 class Header(object):
35 33
36 34 def __init__(self):
37 35 raise NotImplementedError
38 36
39
40 37 def read(self):
41 38
42 39 raise NotImplementedError
@@ -68,13 +65,18 class Header(object):
68 65
69 66 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
70 67 ('Hname','a32'), #Original file name
71 ('Htime',numpy.str_,32), #Date and time when the file was created
72 ('Hoper',numpy.str_,64), #Name of operator who created the file
73 ('Hplace',numpy.str_,128), #Place where the measurements was carried out
74 ('Hdescr',numpy.str_,256), #Description of measurements
68 # Date and time when the file was created
69 ('Htime', numpy.str_, 32),
70 # Name of operator who created the file
71 ('Hoper', numpy.str_, 64),
72 # Place where the measurements was carried out
73 ('Hplace', numpy.str_, 128),
74 # Description of measurements
75 ('Hdescr', numpy.str_, 256),
75 76 ('Hdummy',numpy.str_,512), #Reserved space
76 77 #Main chunk 8bytes
77 ('Msign',numpy.str_,4), #Main chunk signature FZKF or NUIG
78 # Main chunk signature FZKF or NUIG
79 ('Msign', numpy.str_, 4),
78 80 ('MsizeData','<i4'), #Size of data block main chunk
79 81 #Processing DSP parameters 36bytes
80 82 ('PPARsign',numpy.str_,4), #PPAR signature
@@ -82,19 +84,28 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
82 84 ('PPARprf','<i4'), #Pulse repetition frequency
83 85 ('PPARpdr','<i4'), #Pulse duration
84 86 ('PPARsft','<i4'), #FFT length
85 ('PPARavc','<i4'), #Number of spectral (in-coherent) averages
86 ('PPARihp','<i4'), #Number of lowest range gate for moment estimation
87 ('PPARchg','<i4'), #Count for gates for moment estimation
88 ('PPARpol','<i4'), #switch on/off polarimetric measurements. Should be 1.
87 # Number of spectral (in-coherent) averages
88 ('PPARavc', '<i4'),
89 # Number of lowest range gate for moment estimation
90 ('PPARihp', '<i4'),
91 # Count for gates for moment estimation
92 ('PPARchg', '<i4'),
93 # switch on/off polarimetric measurements. Should be 1.
94 ('PPARpol', '<i4'),
89 95 #Service DSP parameters 112bytes
90 ('SPARatt','<i4'), #STC attenuation on the lowest ranges on/off
96 # STC attenuation on the lowest ranges on/off
97 ('SPARatt', '<i4'),
91 98 ('SPARtx','<i4'), #OBSOLETE
92 99 ('SPARaddGain0','<f4'), #OBSOLETE
93 100 ('SPARaddGain1','<f4'), #OBSOLETE
94 ('SPARwnd','<i4'), #Debug only. It normal mode it is 0.
95 ('SPARpos','<i4'), #Delay between sync pulse and tx pulse for phase corr, ns
96 ('SPARadd','<i4'), #"add to pulse" to compensate for delay between the leading edge of driver pulse and envelope of the RF signal.
97 ('SPARlen','<i4'), #Time for measuring txn pulse phase. OBSOLETE
101 # Debug only. It normal mode it is 0.
102 ('SPARwnd', '<i4'),
103 # Delay between sync pulse and tx pulse for phase corr, ns
104 ('SPARpos', '<i4'),
105 # "add to pulse" to compensate for delay between the leading edge of driver pulse and envelope of the RF signal.
106 ('SPARadd', '<i4'),
107 # Time for measuring txn pulse phase. OBSOLETE
108 ('SPARlen', '<i4'),
98 109 ('SPARcal','<i4'), #OBSOLETE
99 110 ('SPARnos','<i4'), #OBSOLETE
100 111 ('SPARof0','<i4'), #detection threshold
@@ -105,19 +116,23 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
105 116 ('SPARtst','<i4'), #OBSOLETE
106 117 ('SPARcor','<i4'), #OBSOLETE
107 118 ('SPARofs','<i4'), #OBSOLETE
108 ('SPARhsn','<i4'), #Hildebrand div noise detection on noise gate
109 ('SPARhsa','<f4'), #Hildebrand div noise detection on all gates
119 # Hildebrand div noise detection on noise gate
120 ('SPARhsn', '<i4'),
121 # Hildebrand div noise detection on all gates
122 ('SPARhsa', '<f4'),
110 123 ('SPARcalibPow_M','<f4'), #OBSOLETE
111 124 ('SPARcalibSNR_M','<f4'), #OBSOLETE
112 125 ('SPARcalibPow_S','<f4'), #OBSOLETE
113 126 ('SPARcalibSNR_S','<f4'), #OBSOLETE
114 ('SPARrawGate1','<i4'), #Lowest range gate for spectra saving Raw_Gate1 >=5
115 ('SPARrawGate2','<i4'), #Number of range gates with atmospheric signal
116 ('SPARraw','<i4'), #flag - IQ or spectra saving on/off
127 # Lowest range gate for spectra saving Raw_Gate1 >=5
128 ('SPARrawGate1', '<i4'),
129 # Number of range gates with atmospheric signal
130 ('SPARrawGate2', '<i4'),
131 # flag - IQ or spectra saving on/off
132 ('SPARraw', '<i4'),
117 133 ('SPARprc','<i4'),]) #flag - Moment estimation switched on/off
118 134
119 135
120
121 136 class FileHeaderMIRA35c(Header):
122 137
123 138 def __init__(self):
@@ -195,7 +210,6 class FileHeaderMIRA35c(Header):
195 210
196 211 '''
197 212
198
199 213 self.Hname= str(header['Hname'][0])
200 214 self.Htime= str(header['Htime'][0])
201 215 self.Hoper= str(header['Hoper'][0])
@@ -272,7 +286,6 class FileHeaderMIRA35c(Header):
272 286 self.Hdescr,
273 287 self.Hdummy)
274 288
275
276 289 header = numpy.array(headerTuple, FILE_HEADER)
277 290 # numpy.array(object, dtype=None, copy=True, order=None, subok=False, ndmin=0)
278 291 header.tofile(fp)
@@ -293,11 +306,13 class FileHeaderMIRA35c(Header):
293 306
294 307 return 1
295 308
309
296 310 SRVI_HEADER = numpy.dtype([
297 ('SignatureSRVI1',numpy.str_,4),#
298 ('SizeOfDataBlock1','<i4'),#
299 ('DataBlockTitleSRVI1',numpy.str_,4),#
300 ('SizeOfSRVI1','<i4'),])#
311 ('SignatureSRVI1', numpy.str_, 4),
312 ('SizeOfDataBlock1', '<i4'),
313 ('DataBlockTitleSRVI1', numpy.str_, 4),
314 ('SizeOfSRVI1', '<i4'), ])
315
301 316
302 317 class SRVIHeader(Header):
303 318 def __init__(self, SignatureSRVI1=0, SizeOfDataBlock1=0, DataBlockTitleSRVI1=0, SizeOfSRVI1=0):
@@ -322,7 +337,7 class SRVIHeader(Header):
322 337
323 338
324 339 SRVI_STRUCTURE = numpy.dtype([
325 ('frame_cnt','<u4'),#
340 ('frame_cnt', '<u4'),
326 341 ('time_t','<u4'), #
327 342 ('tpow','<f4'), #
328 343 ('npw1','<f4'), #
@@ -338,24 +353,20 SRVI_STRUCTURE = numpy.dtype([
338 353 ('azivel','<f4'), #
339 354 ('elvpos','<f4'), #
340 355 ('elvvel','<f4'), #
341 ('northAngle','<f4'), #
356 ('northAngle', '<f4'),
342 357 ('microsec','<u4'), #
343 358 ('azisetvel','<f4'), #
344 359 ('elvsetpos','<f4'), #
345 360 ('RadarConst','<f4'),]) #
346 361
347 362
348
349
350 363 class RecordHeader(Header):
351 364
352
353 365 def __init__(self, frame_cnt=0, time_t= 0, tpow=0, npw1=0, npw2=0,
354 366 cpw1=0, pcw2=0, ps_err=0, te_err=0, rc_err=0, grs1=0,
355 367 grs2=0, azipos=0, azivel=0, elvpos=0, elvvel=0, northangle=0,
356 368 microsec=0, azisetvel=0, elvsetpos=0, RadarConst=0 , RecCounter=0, Off2StartNxtRec=0):
357 369
358
359 370 self.frame_cnt = frame_cnt
360 371 self.dwell = time_t
361 372 self.tpow = tpow
@@ -392,7 +403,7 class RecordHeader(Header):
392 403
393 404 header = numpy.fromfile(fp,SRVI_STRUCTURE,1)
394 405
395 self.frame_cnt = header['frame_cnt'][0]#
406 self.frame_cnt = header['frame_cnt'][0]
396 407 self.time_t = header['time_t'][0] #
397 408 self.tpow = header['tpow'][0] #
398 409 self.npw1 = header['npw1'][0] #
@@ -428,9 +439,9 class RecordHeader(Header):
428 439
429 440 print '=============================================='
430 441
431
432 442 return 1
433 443
444
434 445 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
435 446
436 447 path = None
@@ -441,7 +452,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
441 452 walk = None
442 453 isConfig = False
443 454
444
445 455 fileList= None
446 456
447 457 #metadata
@@ -453,8 +463,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
453 463 data= None
454 464 utctime= None
455 465
456
457
458 466 def __init__(self, **kwargs):
459 467
460 468 #Eliminar de la base la herencia
@@ -498,7 +506,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
498 506 self.dataOut.COFA = []
499 507 self.dataOut.noise = 0
500 508
501
502 509 def Files2Read(self, fp):
503 510 '''
504 511 Function that indicates the number of .fdt files that exist in the folder to be read.
@@ -506,8 +513,10 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
506 513 '''
507 514 #self.__checkPath()
508 515
509 ListaData=os.listdir(fp) #Gets the list of files within the fp address
510 ListaData=sorted(ListaData) #Sort the list of files from least to largest by names
516 # Gets the list of files within the fp address
517 ListaData = os.listdir(fp)
518 # Sort the list of files from least to largest by names
519 ListaData = sorted(ListaData)
511 520 nFiles=0 #File Counter
512 521 FileList=[] #A list is created that will contain the .fdt files
513 522 for IndexFile in ListaData :
@@ -520,7 +529,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
520 529
521 530 self.filenameList=FileList #List of files from least to largest by names
522 531
523
524 532 def run(self, **kwargs):
525 533 '''
526 534 This method will be the one that will initiate the data entry, will be called constantly.
@@ -533,7 +541,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
533 541
534 542 self.getData()
535 543
536
537 544 def setup(self, path=None,
538 545 startDate=None,
539 546 endDate=None,
@@ -557,7 +564,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
557 564
558 565 pass
559 566
560
561 567 def getData(self):
562 568 '''
563 569 Before starting this function, you should check that there is still an unread file,
@@ -584,10 +590,8 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
584 590 self.dataOut.data_spc = self.dataOut.data_spc+self.dataOut.noise
585 591 #print 'self.dataOut.noise',self.dataOut.noise
586 592
587
588 593 return self.dataOut.data_spc
589 594
590
591 595 def readFile(self,fp):
592 596 '''
593 597 You must indicate if you are reading in Online or Offline mode and load the
@@ -604,7 +608,8 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
604 608
605 609 if self.fileSelector < len(self.filenameList):
606 610
607 self.fpFile=str(fp)+'/'+str(self.filenameList[self.fileSelector])
611 self.fpFile = str(fp) + '/' + \
612 str(self.filenameList[self.fileSelector])
608 613
609 614 if self.nextfileflag==True:
610 615 self.fp = open(self.fpFile,"rb")
@@ -612,18 +617,15 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
612 617
613 618 '''HERE STARTING THE FILE READING'''
614 619
615
616 620 self.fheader = FileHeaderMIRA35c()
617 621 self.fheader.FHread(self.fp) #Bltr FileHeader Reading
618 622
619
620 623 self.SPARrawGate1 = self.fheader.SPARrawGate1
621 624 self.SPARrawGate2 = self.fheader.SPARrawGate2
622 625 self.Num_Hei = self.SPARrawGate2 - self.SPARrawGate1
623 626 self.Num_Bins = self.fheader.PPARsft
624 627 self.dataOut.nFFTPoints = self.fheader.PPARsft
625 628
626
627 629 self.Num_inCoh = self.fheader.PPARavc
628 630 self.dataOut.PRF = self.fheader.PPARprf
629 631 self.dataOut.frequency = 34.85*10**9
@@ -633,14 +635,12 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
633 635 pulse_width = self.fheader.PPARpdr * 10**-9
634 636 self.__deltaHeigth = 0.5 * SPEED_OF_LIGHT * pulse_width
635 637
636 self.data_spc = numpy.zeros((self.Num_Hei, self.Num_Bins,2))#
638 self.data_spc = numpy.zeros((self.Num_Hei, self.Num_Bins, 2))
637 639 self.dataOut.HSDV = numpy.zeros((self.Num_Hei, 2))
638 640
639 641 self.Ze = numpy.zeros(self.Num_Hei)
640 642 self.ETA = numpy.zeros(([2,self.Num_Hei]))
641 643
642
643
644 644 self.readBlock() #Block reading
645 645
646 646 else:
@@ -650,8 +650,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
650 650 self.FileHeaderFlag == True
651 651 return 0
652 652
653
654
655 653 def readBlock(self):
656 654 '''
657 655 It should be checked if the block has data, if it is not passed to the next file.
@@ -670,8 +668,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
670 668 else :
671 669 self.FirstPoint = 1180
672 670
673
674
675 671 self.srviHeader = SRVIHeader()
676 672
677 673 self.srviHeader.SRVIread(self.fp) #Se obtiene la cabecera del SRVI
@@ -682,7 +678,8 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
682 678 print 'blocksize == 148 bug'
683 679 jump = numpy.fromfile(self.fp,[('jump',numpy.str_,140)] ,1)
684 680
685 self.srviHeader.SRVIread(self.fp) #Se obtiene la cabecera del SRVI
681 # Se obtiene la cabecera del SRVI
682 self.srviHeader.SRVIread(self.fp)
686 683
687 684 if not self.srviHeader.SizeOfSRVI1:
688 685 self.fileSelector+=1
@@ -696,7 +693,6 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
696 693 npw1 = self.recordheader.npw1
697 694 npw2 = self.recordheader.npw2
698 695
699
700 696 self.dataOut.channelList = range(1)
701 697 self.dataOut.nIncohInt = self.Num_inCoh
702 698 self.dataOut.nProfiles = self.Num_Bins
@@ -706,53 +702,65 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
706 702 self.dataOut.timeZone=0
707 703
708 704 self.dataOut.outputInterval = self.dataOut.getTimeInterval()
709 self.dataOut.heightList = self.SPARrawGate1*self.__deltaHeigth + numpy.array(range(self.Num_Hei))*self.__deltaHeigth
710
711
705 self.dataOut.heightList = self.SPARrawGate1 * self.__deltaHeigth + \
706 numpy.array(range(self.Num_Hei)) * self.__deltaHeigth
712 707
713 708 self.HSDVsign = numpy.fromfile( self.fp, [('HSDV',numpy.str_,4)],1)
714 709 self.SizeHSDV = numpy.fromfile( self.fp, [('SizeHSDV','<i4')],1)
715 self.HSDV_Co = numpy.fromfile( self.fp, [('HSDV_Co','<f4')],self.Num_Hei)
716 self.HSDV_Cx = numpy.fromfile( self.fp, [('HSDV_Cx','<f4')],self.Num_Hei)
710 self.HSDV_Co = numpy.fromfile(
711 self.fp, [('HSDV_Co', '<f4')], self.Num_Hei)
712 self.HSDV_Cx = numpy.fromfile(
713 self.fp, [('HSDV_Cx', '<f4')], self.Num_Hei)
717 714
718 715 self.COFAsign = numpy.fromfile( self.fp, [('COFA',numpy.str_,4)],1)
719 716 self.SizeCOFA = numpy.fromfile( self.fp, [('SizeCOFA','<i4')],1)
720 self.COFA_Co = numpy.fromfile( self.fp, [('COFA_Co','<f4')],self.Num_Hei)
721 self.COFA_Cx = numpy.fromfile( self.fp, [('COFA_Cx','<f4')],self.Num_Hei)
717 self.COFA_Co = numpy.fromfile(
718 self.fp, [('COFA_Co', '<f4')], self.Num_Hei)
719 self.COFA_Cx = numpy.fromfile(
720 self.fp, [('COFA_Cx', '<f4')], self.Num_Hei)
722 721
723 self.ZSPCsign = numpy.fromfile(self.fp, [('ZSPCsign',numpy.str_,4)],1)
722 self.ZSPCsign = numpy.fromfile(
723 self.fp, [('ZSPCsign', numpy.str_, 4)], 1)
724 724 self.SizeZSPC = numpy.fromfile(self.fp, [('SizeZSPC','<i4')],1)
725 725
726 726 self.dataOut.HSDV[0]=self.HSDV_Co[:][0]
727 727 self.dataOut.HSDV[1]=self.HSDV_Cx[:][0]
728 728
729 729 for irg in range(self.Num_Hei):
730 nspc = numpy.fromfile(self.fp, [('nspc','int16')],1)[0][0] # Number of spectral sub pieces containing significant power
730 # Number of spectral sub pieces containing significant power
731 nspc = numpy.fromfile(self.fp, [('nspc', 'int16')], 1)[0][0]
731 732
732 733 for k in range(nspc):
733 binIndex = numpy.fromfile(self.fp, [('binIndex','int16')],1)[0][0] # Index of the spectral bin where the piece is beginning
734 nbins = numpy.fromfile(self.fp, [('nbins','int16')],1)[0][0] # Number of bins of the piece
734 # Index of the spectral bin where the piece is beginning
735 binIndex = numpy.fromfile(
736 self.fp, [('binIndex', 'int16')], 1)[0][0]
737 nbins = numpy.fromfile(self.fp, [('nbins', 'int16')], 1)[
738 0][0] # Number of bins of the piece
735 739
736 740 #Co_Channel
737 jbin = numpy.fromfile(self.fp, [('jbin','uint16')],nbins)[0][0] # Spectrum piece to be normaliced
738 jmax = numpy.fromfile(self.fp, [('jmax','float32')],1)[0][0] # Maximun piece to be normaliced
739
741 jbin = numpy.fromfile(self.fp, [('jbin', 'uint16')], nbins)[
742 0][0] # Spectrum piece to be normaliced
743 jmax = numpy.fromfile(self.fp, [('jmax', 'float32')], 1)[
744 0][0] # Maximun piece to be normaliced
740 745
741 self.data_spc[irg,binIndex:binIndex+nbins,0] = self.data_spc[irg,binIndex:binIndex+nbins,0]+jbin/65530.*jmax
746 self.data_spc[irg, binIndex:binIndex + nbins, 0] = self.data_spc[irg,
747 binIndex:binIndex + nbins, 0] + jbin / 65530. * jmax
742 748
743 749 #Cx_Channel
744 jbin = numpy.fromfile(self.fp, [('jbin','uint16')],nbins)[0][0]
750 jbin = numpy.fromfile(
751 self.fp, [('jbin', 'uint16')], nbins)[0][0]
745 752 jmax = numpy.fromfile(self.fp, [('jmax','float32')],1)[0][0]
746 753
747
748 self.data_spc[irg,binIndex:binIndex+nbins,1] = self.data_spc[irg,binIndex:binIndex+nbins,1]+jbin/65530.*jmax
754 self.data_spc[irg, binIndex:binIndex + nbins, 1] = self.data_spc[irg,
755 binIndex:binIndex + nbins, 1] + jbin / 65530. * jmax
749 756
750 757 for bin in range(self.Num_Bins):
751 758
752 self.data_spc[:,bin,0] = self.data_spc[:,bin,0] - self.dataOut.HSDV[:,0]
753
754 self.data_spc[:,bin,1] = self.data_spc[:,bin,1] - self.dataOut.HSDV[:,1]
759 self.data_spc[:, bin, 0] = self.data_spc[:,
760 bin, 0] - self.dataOut.HSDV[:, 0]
755 761
762 self.data_spc[:, bin, 1] = self.data_spc[:,
763 bin, 1] - self.dataOut.HSDV[:, 1]
756 764
757 765 numpy.set_printoptions(threshold='nan')
758 766
@@ -787,17 +795,8 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
787 795 #shift = self.Num_Bins/2 + fix
788 796 #self.data_spc = numpy.array([ self.data_spc[: , self.Num_Bins-shift+1: , :] , self.data_spc[: , 0:self.Num_Bins-shift , :]])
789 797
790
791
792 798 '''Block Reading, the Block Data is received and Reshape is used to give it
793 799 shape.
794 800 '''
795 801
796 802 self.PointerReader = self.fp.tell()
797
798
799
800
801
802
803 No newline at end of file
@@ -17,7 +17,6 from functools import partial
17 17 import time
18 18 #from sklearn.cluster import KMeans
19 19
20 import matplotlib.pyplot as plt
21 20
22 21 from scipy.optimize import fmin_l_bfgs_b #optimize with bounds on state papameters
23 22 from jroproc_base import ProcessingUnit, Operation
General Comments 0
You need to be logged in to leave comments. Login now