##// END OF EJS Templates
schain
RSS
sin matplotlib de los modulos de claire
José Chávez -
r1075:eb67833bff79
parent child
Show More
Side by Side Unified
Context file:
r1075:eb67833bff79 Loading diff...:
Show file before | Show file after | Hide comments
@@ -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 from utils import * No newline at end of file
12 from utils import *
Show file before | Show file after | Hide comments
@@ -4,32 +4,36 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')
11 #Qt4Agg', 'GTK', 'GTKAgg', 'ps', 'agg', 'cairo', 'MacOSX', 'GTKCairo', 'WXAgg', 'template', 'TkAgg', 'GTK3Cairo', 'GTK3Agg', 'svg', 'WebAgg', 'CocoaAgg', 'emf', 'gdk', 'WX'
11 # Qt4Agg', 'GTK', 'GTKAgg', 'ps', 'agg', 'cairo', 'MacOSX', 'GTKCairo', 'WXAgg', 'template', 'TkAgg', 'GTK3Cairo', 'GTK3Agg', 'svg', 'WebAgg', 'CocoaAgg', 'emf', 'gdk', 'WX'
12 12 import matplotlib.pyplot
13 13
14 14 from mpl_toolkits.axes_grid1 import make_axes_locatable
15 15 from matplotlib.ticker import FuncFormatter, LinearLocator
16 16
17 17 ###########################################
18 #Actualizacion de las funciones del driver
18 # Actualizacion de las funciones del driver
19 19 ###########################################
20 20
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
28 def createFigure(id, wintitle, width, height, facecolor="w", show=True, dpi = 80):
30
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,10 +43,11 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 # matplotlib.pyplot.ioff()
45 # matplotlib.pyplot.pause(0)
49 # matplotlib.pyplot.ioff()
50 # matplotlib.pyplot.pause(0)
46 51
47 52 if show:
48 53 matplotlib.pyplot.show()
@@ -60,45 +65,52 def closeFigure(show=False, fig=None):
60 65
61 66 return
62 67
68
63 69 def saveFigure(fig, filename):
64 70
65 # matplotlib.pyplot.ioff()
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()
84 94 matplotlib.pyplot.figure(fig.number)
85 95 axes = matplotlib.pyplot.subplot2grid((nrow, ncol),
86 (xpos, ypos),
87 colspan=colspan,
88 rowspan=rowspan,
89 polar=polar)
96 (xpos, ypos),
97 colspan=colspan,
98 rowspan=rowspan,
99 polar=polar)
90 100
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,
97 xy = (.1, .99),
98 xycoords = 'figure fraction',
99 horizontalalignment = 'left',
100 verticalalignment = 'top',
101 fontsize = 10)
108 xy=(.1, .99),
109 xycoords='figure fraction',
110 horizontalalignment='left',
111 verticalalignment='top',
112 fontsize=10)
113
102 114
103 115 def printLabels(ax, xlabel, ylabel, title):
104 116
@@ -106,11 +118,11 def printLabels(ax, xlabel, ylabel, title):
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 grid=None,color='blue'):
113
125 grid=None, color='blue'):
114 126 """
115 127
116 128 Input:
@@ -119,18 +131,19 def createPline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title=''
119 131
120 132 matplotlib.pyplot.ioff()
121 133
122 ax.set_xlim([xmin,xmax])
123 ax.set_ylim([ymin,ymax])
134 ax.set_xlim([xmin, xmax])
135 ax.set_ylim([ymin, ymax])
124 136
125 137 printLabels(ax, xlabel, ylabel, title)
126 138
127 139 ######################################################
128 if (xmax-xmin)<=1:
129 xtickspos = numpy.linspace(xmin,xmax,nxticks)
130 xtickspos = numpy.array([float("%.1f"%i) for i in xtickspos])
140 if (xmax - xmin) <= 1:
141 xtickspos = numpy.linspace(xmin, xmax, nxticks)
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,9 +181,11 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 ax.lines[idline].set_data(x,y)
187 ax.lines[idline].set_data(x, y)
188
174 189
175 190 def pline(iplot, x, y, xlabel='', ylabel='', title=''):
176 191
@@ -180,14 +195,15 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 ax.plot(x,y,color=color,linestyle=linestyle,lw=lw)
201 ax.plot(x, y, color=color, linestyle=linestyle, lw=lw)
186 202
187 203
188 204 def createPcolor(ax, x, y, z, xmin, xmax, ymin, ymax, zmin, zmax,
189 xlabel='', ylabel='', title='', ticksize = 9,
190 colormap='jet',cblabel='', cbsize="5%",
205 xlabel='', ylabel='', title='', ticksize=9,
206 colormap='jet', cblabel='', cbsize="5%",
191 207 XAxisAsTime=False):
192 208
193 209 matplotlib.pyplot.ioff()
@@ -197,16 +213,16 def createPcolor(ax, x, y, z, xmin, xmax, ymin, ymax, zmin, zmax,
197 213 fig = ax.get_figure()
198 214 fig.add_axes(ax_cb)
199 215
200 ax.set_xlim([xmin,xmax])
201 ax.set_ylim([ymin,ymax])
216 ax.set_xlim([xmin, xmax])
217 ax.set_ylim([ymin, ymax])
202 218
203 219 printLabels(ax, xlabel, ylabel, title)
204 220
205 221 z = numpy.ma.masked_invalid(z)
206 cmap=matplotlib.pyplot.get_cmap(colormap)
222 cmap = matplotlib.pyplot.get_cmap(colormap)
207 223 cmap.set_bad('black', 1.)
208 imesh = ax.pcolormesh(x,y,z.T, vmin=zmin, vmax=zmax, cmap=cmap)
209 cb = matplotlib.pyplot.colorbar(imesh, cax=ax_cb)
224 imesh = ax.pcolormesh(x, y, z.T, vmin=zmin, vmax=zmax, cmap=cmap)
225 cb = matplotlib.pyplot.colorbar(imesh, cax=ax_cb)
210 226 cb.set_label(cblabel)
211 227
212 228 # for tl in ax_cb.get_yticklabels():
@@ -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
@@ -263,17 +284,16 def addpcolorbuffer(ax, x, y, z, zmin, zmax, xlabel='', ylabel='', title='', col
263 284
264 285 z = numpy.ma.masked_invalid(z)
265 286
266 cmap=matplotlib.pyplot.get_cmap(colormap)
287 cmap = matplotlib.pyplot.get_cmap(colormap)
267 288 cmap.set_bad('black', 1.)
268 289
290 ax.pcolormesh(x, y, z.T, vmin=zmin, vmax=zmax, cmap=cmap)
269 291
270 ax.pcolormesh(x,y,z.T,vmin=zmin,vmax=zmax, cmap=cmap)
271 292
272 293 def createPmultiline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
273 ticksize=9, xtick_visible=True, ytick_visible=True,
274 nxticks=4, nyticks=10,
275 grid=None):
276
294 ticksize=9, xtick_visible=True, ytick_visible=True,
295 nxticks=4, nyticks=10,
296 grid=None):
277 297 """
278 298
279 299 Input:
@@ -285,11 +305,12 def createPmultiline(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', tit
285 305 lines = ax.plot(x.T, y)
286 306 leg = ax.legend(lines, legendlabels, loc='upper right')
287 307 leg.get_frame().set_alpha(0.5)
288 ax.set_xlim([xmin,xmax])
289 ax.set_ylim([ymin,ymax])
308 ax.set_xlim([xmin, xmax])
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():
@@ -332,13 +353,13 def pmultiline(iplot, x, y, xlabel='', ylabel='', title=''):
332 353
333 354 for i in range(len(ax.lines)):
334 355 line = ax.lines[i]
335 line.set_data(x[i,:],y)
356 line.set_data(x[i, :], y)
336 357
337 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
338 ticksize=9, xtick_visible=True, ytick_visible=True,
339 nxticks=4, nyticks=10, marker='.', markersize=10, linestyle="None",
340 grid=None, XAxisAsTime=False):
341 358
359 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='', title='', legendlabels=None,
360 ticksize=9, xtick_visible=True, ytick_visible=True,
361 nxticks=4, nyticks=10, marker='.', markersize=10, linestyle="None",
362 grid=None, XAxisAsTime=False):
342 363 """
343 364
344 365 Input:
@@ -355,10 +376,11 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 ax.set_xlim([xmin,xmax])
361 ax.set_ylim([ymin,ymax])
382 ax.set_xlim([xmin, xmax])
383 ax.set_ylim([ymin, ymax])
362 384 printLabels(ax, xlabel, ylabel, title)
363 385
364 386 # xtickspos = numpy.arange(nxticks)*int((xmax-xmin)/(nxticks)) + int(xmin)
@@ -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
@@ -409,19 +433,20 def pmultilineyaxis(iplot, x, y, xlabel='', ylabel='', title=''):
409 433
410 434 for i in range(len(ax.lines)):
411 435 line = ax.lines[i]
412 line.set_data(x,y[i,:])
436 line.set_data(x, y[i, :])
437
413 438
414 439 def createPolar(ax, x, y,
415 xlabel='', ylabel='', title='', ticksize = 9,
416 colormap='jet',cblabel='', cbsize="5%",
417 XAxisAsTime=False):
440 xlabel='', ylabel='', title='', ticksize=9,
441 colormap='jet', cblabel='', cbsize="5%",
442 XAxisAsTime=False):
418 443
419 444 matplotlib.pyplot.ioff()
420 445
421 ax.plot(x,y,'bo', markersize=5)
446 ax.plot(x, y, 'bo', markersize=5)
422 447 # ax.set_rmax(90)
423 ax.set_ylim(0,90)
424 ax.set_yticks(numpy.arange(0,90,20))
448 ax.set_ylim(0, 90)
449 ax.set_yticks(numpy.arange(0, 90, 20))
425 450 # ax.text(0, -110, ylabel, rotation='vertical', va ='center', ha = 'center' ,size='11')
426 451 # ax.text(0, 50, ylabel, rotation='vertical', va ='center', ha = 'left' ,size='11')
427 452 # ax.text(100, 100, 'example', ha='left', va='center', rotation='vertical')
@@ -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)
Show file before | Show file after | Hide comments
@@ -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,11 +7,9 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 from scipy import asarray as ar,exp
12 from scipy import asarray as ar, exp
14 13 from scipy import stats
15 14
16 15 from duplicity.path import Path
@@ -31,113 +30,130 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")
35
36
37 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
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
43 ('Hdummy',numpy.str_,512), #Reserved space
44 #Main chunk
45 ('Msign','<i4'), #Main chunk signature FZKF or NUIG
46 ('MsizeData','<i4'), #Size of data block main chunk
47 #Processing DSP parameters
48 ('PPARsign','<i4'), #PPAR signature
49 ('PPARsize','<i4'), #PPAR size of block
50 ('PPARprf','<i4'), #Pulse repetition frequency
51 ('PPARpdr','<i4'), #Pulse duration
52 ('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 #Service DSP parameters
58 ('SPARatt','<i4'), #STC attenuation on the lowest ranges on/off
59 ('SPARtx','<i4'), #OBSOLETE
60 ('SPARaddGain0','<f4'), #OBSOLETE
61 ('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
66 ('SPARcal','<i4'), #OBSOLETE
67 ('SPARnos','<i4'), #OBSOLETE
68 ('SPARof0','<i4'), #detection threshold
69 ('SPARof1','<i4'), #OBSOLETE
70 ('SPARswt','<i4'), #2nd moment estimation threshold
71 ('SPARsum','<i4'), #OBSOLETE
72 ('SPARosc','<i4'), #flag Oscillosgram mode
73 ('SPARtst','<i4'), #OBSOLETE
74 ('SPARcor','<i4'), #OBSOLETE
75 ('SPARofs','<i4'), #OBSOLETE
76 ('SPARhsn','<i4'), #Hildebrand div noise detection on noise gate
77 ('SPARhsa','<f4'), #Hildebrand div noise detection on all gates
78 ('SPARcalibPow_M','<f4'), #OBSOLETE
79 ('SPARcalibSNR_M','<f4'), #OBSOLETE
80 ('SPARcalibPow_S','<f4'), #OBSOLETE
81 ('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
85 ('SPARprc','<i4'),]) #flag - Moment estimation switched on/off
86
87
88
89 self.Hname= None
90 self.Htime= None
91 self.Hoper= None
92 self.Hplace= None
93 self.Hdescr= None
94 self.Hdummy= None
95
96 self.Msign=None
97 self.MsizeData=None
98
99 self.PPARsign=None
100 self.PPARsize=None
101 self.PPARprf=None
102 self.PPARpdr=None
103 self.PPARsft=None
104 self.PPARavc=None
105 self.PPARihp=None
106 self.PPARchg=None
107 self.PPARpol=None
108 #Service DSP parameters
109 self.SPARatt=None
110 self.SPARtx=None
111 self.SPARaddGain0=None
112 self.SPARaddGain1=None
113 self.SPARwnd=None
114 self.SPARpos=None
115 self.SPARadd=None
116 self.SPARlen=None
117 self.SPARcal=None
118 self.SPARnos=None
119 self.SPARof0=None
120 self.SPARof1=None
121 self.SPARswt=None
122 self.SPARsum=None
123 self.SPARosc=None
124 self.SPARtst=None
125 self.SPARcor=None
126 self.SPARofs=None
127 self.SPARhsn=None
128 self.SPARhsa=None
129 self.SPARcalibPow_M=None
130 self.SPARcalibSNR_M=None
131 self.SPARcalibPow_S=None
132 self.SPARcalibSNR_S=None
133 self.SPARrawGate1=None
134 self.SPARrawGate2=None
135 self.SPARraw=None
136 self.SPARprc=None
137
138
139
140 header = numpy.fromfile(fp, FILE_HEADER,1)
33 startFp = open(
34 '/home/erick/Documents/MIRA35C/20160117/20160117_0000.zspc', "rb")
35
36
37 FILE_HEADER = numpy.dtype([ # HEADER 1024bytes
38 ('Hname', numpy.str_, 32), # Original file name
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),
47 ('Hdummy', numpy.str_, 512), # Reserved space
48 # Main chunk
49 ('Msign', '<i4'), # Main chunk signature FZKF or NUIG
50 ('MsizeData', '<i4'), # Size of data block main chunk
51 # Processing DSP parameters
52 ('PPARsign', '<i4'), # PPAR signature
53 ('PPARsize', '<i4'), # PPAR size of block
54 ('PPARprf', '<i4'), # Pulse repetition frequency
55 ('PPARpdr', '<i4'), # Pulse duration
56 ('PPARsft', '<i4'), # FFT length
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'),
65 # Service DSP parameters
66 # STC attenuation on the lowest ranges on/off
67 ('SPARatt', '<i4'),
68 ('SPARtx', '<i4'), # OBSOLETE
69 ('SPARaddGain0', '<f4'), # OBSOLETE
70 ('SPARaddGain1', '<f4'), # 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'),
79 ('SPARcal', '<i4'), # OBSOLETE
80 ('SPARnos', '<i4'), # OBSOLETE
81 ('SPARof0', '<i4'), # detection threshold
82 ('SPARof1', '<i4'), # OBSOLETE
83 ('SPARswt', '<i4'), # 2nd moment estimation threshold
84 ('SPARsum', '<i4'), # OBSOLETE
85 ('SPARosc', '<i4'), # flag Oscillosgram mode
86 ('SPARtst', '<i4'), # OBSOLETE
87 ('SPARcor', '<i4'), # OBSOLETE
88 ('SPARofs', '<i4'), # OBSOLETE
89 # Hildebrand div noise detection on noise gate
90 ('SPARhsn', '<i4'),
91 # Hildebrand div noise detection on all gates
92 ('SPARhsa', '<f4'),
93 ('SPARcalibPow_M', '<f4'), # OBSOLETE
94 ('SPARcalibSNR_M', '<f4'), # OBSOLETE
95 ('SPARcalibPow_S', '<f4'), # OBSOLETE
96 ('SPARcalibSNR_S', '<f4'), # OBSOLETE
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'),
103 ('SPARprc', '<i4'), ]) # flag - Moment estimation switched on/off
104
105
106 self.Hname = None
107 self.Htime = None
108 self.Hoper = None
109 self.Hplace = None
110 self.Hdescr = None
111 self.Hdummy = None
112
113 self.Msign = None
114 self.MsizeData = None
115
116 self.PPARsign = None
117 self.PPARsize = None
118 self.PPARprf = None
119 self.PPARpdr = None
120 self.PPARsft = None
121 self.PPARavc = None
122 self.PPARihp = None
123 self.PPARchg = None
124 self.PPARpol = None
125 # Service DSP parameters
126 self.SPARatt = None
127 self.SPARtx = None
128 self.SPARaddGain0 = None
129 self.SPARaddGain1 = None
130 self.SPARwnd = None
131 self.SPARpos = None
132 self.SPARadd = None
133 self.SPARlen = None
134 self.SPARcal = None
135 self.SPARnos = None
136 self.SPARof0 = None
137 self.SPARof1 = None
138 self.SPARswt = None
139 self.SPARsum = None
140 self.SPARosc = None
141 self.SPARtst = None
142 self.SPARcor = None
143 self.SPARofs = None
144 self.SPARhsn = None
145 self.SPARhsa = None
146 self.SPARcalibPow_M = None
147 self.SPARcalibSNR_M = None
148 self.SPARcalibPow_S = None
149 self.SPARcalibSNR_S = None
150 self.SPARrawGate1 = None
151 self.SPARrawGate2 = None
152 self.SPARraw = None
153 self.SPARprc = None
154
155
156 header = numpy.fromfile(fp, FILE_HEADER, 1)
141 157 ''' numpy.fromfile(file, dtype, count, sep='')
142 158 file : file or str
143 159 Open file object or filename.
@@ -157,95 +173,93 header = numpy.fromfile(fp, FILE_HEADER,1)
157 173
158 174 '''
159 175
160 Hname= str(header['Hname'][0])
161 Htime= str(header['Htime'][0])
162 Hoper= str(header['Hoper'][0])
163 Hplace= str(header['Hplace'][0])
164 Hdescr= str(header['Hdescr'][0])
165 Hdummy= str(header['Hdummy'][0])
166
167 Msign=header['Msign'][0]
168 MsizeData=header['MsizeData'][0]
169
170 PPARsign=header['PPARsign'][0]
171 PPARsize=header['PPARsize'][0]
172 PPARprf=header['PPARprf'][0]
173 PPARpdr=header['PPARpdr'][0]
174 PPARsft=header['PPARsft'][0]
175 PPARavc=header['PPARavc'][0]
176 PPARihp=header['PPARihp'][0]
177 PPARchg=header['PPARchg'][0]
178 PPARpol=header['PPARpol'][0]
179 #Service DSP parameters
180 SPARatt=header['SPARatt'][0]
181 SPARtx=header['SPARtx'][0]
182 SPARaddGain0=header['SPARaddGain0'][0]
183 SPARaddGain1=header['SPARaddGain1'][0]
184 SPARwnd=header['SPARwnd'][0]
185 SPARpos=header['SPARpos'][0]
186 SPARadd=header['SPARadd'][0]
187 SPARlen=header['SPARlen'][0]
188 SPARcal=header['SPARcal'][0]
189 SPARnos=header['SPARnos'][0]
190 SPARof0=header['SPARof0'][0]
191 SPARof1=header['SPARof1'][0]
192 SPARswt=header['SPARswt'][0]
193 SPARsum=header['SPARsum'][0]
194 SPARosc=header['SPARosc'][0]
195 SPARtst=header['SPARtst'][0]
196 SPARcor=header['SPARcor'][0]
197 SPARofs=header['SPARofs'][0]
198 SPARhsn=header['SPARhsn'][0]
199 SPARhsa=header['SPARhsa'][0]
200 SPARcalibPow_M=header['SPARcalibPow_M'][0]
201 SPARcalibSNR_M=header['SPARcalibSNR_M'][0]
202 SPARcalibPow_S=header['SPARcalibPow_S'][0]
203 SPARcalibSNR_S=header['SPARcalibSNR_S'][0]
204 SPARrawGate1=header['SPARrawGate1'][0]
205 SPARrawGate2=header['SPARrawGate2'][0]
206 SPARraw=header['SPARraw'][0]
207 SPARprc=header['SPARprc'][0]
208
209
210
211 SRVI_STRUCTURE = numpy.dtype([
212 ('frame_cnt','<u4'),#
213 ('time_t','<u4'), #
214 ('tpow','<f4'), #
215 ('npw1','<f4'), #
216 ('npw2','<f4'), #
217 ('cpw1','<f4'), #
218 ('pcw2','<f4'), #
219 ('ps_err','<u4'), #
220 ('te_err','<u4'), #
221 ('rc_err','<u4'), #
222 ('grs1','<u4'), #
223 ('grs2','<u4'), #
224 ('azipos','<f4'), #
225 ('azivel','<f4'), #
226 ('elvpos','<f4'), #
227 ('elvvel','<f4'), #
228 ('northAngle','<f4'), #
229 ('microsec','<u4'), #
230 ('azisetvel','<f4'), #
231 ('elvsetpos','<f4'), #
232 ('RadarConst','<f4'),]) #
233
234 JUMP_STRUCTURE = numpy.dtype([
235 ('jump','<u140'),#
236 ('SizeOfDataBlock1',numpy.str_,32),#
237 ('jump','<i4'),#
238 ('DataBlockTitleSRVI1',numpy.str_,32),#
239 ('SizeOfSRVI1','<i4'),])#
240
241
242
243 #frame_cnt=0, time_t= 0, tpow=0, npw1=0, npw2=0,
244 #cpw1=0, pcw2=0, ps_err=0, te_err=0, rc_err=0, grs1=0,
245 #grs2=0, azipos=0, azivel=0, elvpos=0, elvvel=0, northangle=0,
246 #microsec=0, azisetvel=0, elvsetpos=0, RadarConst=0
247
248
176 Hname = str(header['Hname'][0])
177 Htime = str(header['Htime'][0])
178 Hoper = str(header['Hoper'][0])
179 Hplace = str(header['Hplace'][0])
180 Hdescr = str(header['Hdescr'][0])
181 Hdummy = str(header['Hdummy'][0])
182
183 Msign = header['Msign'][0]
184 MsizeData = header['MsizeData'][0]
185
186 PPARsign = header['PPARsign'][0]
187 PPARsize = header['PPARsize'][0]
188 PPARprf = header['PPARprf'][0]
189 PPARpdr = header['PPARpdr'][0]
190 PPARsft = header['PPARsft'][0]
191 PPARavc = header['PPARavc'][0]
192 PPARihp = header['PPARihp'][0]
193 PPARchg = header['PPARchg'][0]
194 PPARpol = header['PPARpol'][0]
195 # Service DSP parameters
196 SPARatt = header['SPARatt'][0]
197 SPARtx = header['SPARtx'][0]
198 SPARaddGain0 = header['SPARaddGain0'][0]
199 SPARaddGain1 = header['SPARaddGain1'][0]
200 SPARwnd = header['SPARwnd'][0]
201 SPARpos = header['SPARpos'][0]
202 SPARadd = header['SPARadd'][0]
203 SPARlen = header['SPARlen'][0]
204 SPARcal = header['SPARcal'][0]
205 SPARnos = header['SPARnos'][0]
206 SPARof0 = header['SPARof0'][0]
207 SPARof1 = header['SPARof1'][0]
208 SPARswt = header['SPARswt'][0]
209 SPARsum = header['SPARsum'][0]
210 SPARosc = header['SPARosc'][0]
211 SPARtst = header['SPARtst'][0]
212 SPARcor = header['SPARcor'][0]
213 SPARofs = header['SPARofs'][0]
214 SPARhsn = header['SPARhsn'][0]
215 SPARhsa = header['SPARhsa'][0]
216 SPARcalibPow_M = header['SPARcalibPow_M'][0]
217 SPARcalibSNR_M = header['SPARcalibSNR_M'][0]
218 SPARcalibPow_S = header['SPARcalibPow_S'][0]
219 SPARcalibSNR_S = header['SPARcalibSNR_S'][0]
220 SPARrawGate1 = header['SPARrawGate1'][0]
221 SPARrawGate2 = header['SPARrawGate2'][0]
222 SPARraw = header['SPARraw'][0]
223 SPARprc = header['SPARprc'][0]
224
225
226 SRVI_STRUCTURE = numpy.dtype([
227 ('frame_cnt', '<u4'),
228 ('time_t', '<u4'), #
229 ('tpow', '<f4'), #
230 ('npw1', '<f4'), #
231 ('npw2', '<f4'), #
232 ('cpw1', '<f4'), #
233 ('pcw2', '<f4'), #
234 ('ps_err', '<u4'), #
235 ('te_err', '<u4'), #
236 ('rc_err', '<u4'), #
237 ('grs1', '<u4'), #
238 ('grs2', '<u4'), #
239 ('azipos', '<f4'), #
240 ('azivel', '<f4'), #
241 ('elvpos', '<f4'), #
242 ('elvvel', '<f4'), #
243 ('northAngle', '<f4'),
244 ('microsec', '<u4'), #
245 ('azisetvel', '<f4'), #
246 ('elvsetpos', '<f4'), #
247 ('RadarConst', '<f4'), ]) #
248
249 JUMP_STRUCTURE = numpy.dtype([
250 ('jump', '<u140'),
251 ('SizeOfDataBlock1', numpy.str_, 32),
252 ('jump', '<i4'),
253 ('DataBlockTitleSRVI1', numpy.str_, 32),
254 ('SizeOfSRVI1', '<i4'), ])
255
256
257 # frame_cnt=0, time_t= 0, tpow=0, npw1=0, npw2=0,
258 # cpw1=0, pcw2=0, ps_err=0, te_err=0, rc_err=0, grs1=0,
259 # grs2=0, azipos=0, azivel=0, elvpos=0, elvvel=0, northangle=0,
260 # microsec=0, azisetvel=0, elvsetpos=0, RadarConst=0
261
262
249 263 frame_cnt = frame_cnt
250 264 dwell = time_t
251 265 tpow = tpow
@@ -267,23 +281,22 microsec = microsec
267 281 azisetvel = azisetvel
268 282 elvsetpos = elvsetpos
269 283 RadarConst5 = RadarConst
270
271 284
272 285
273 #print fp
274 #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 #startFp = open(fp,"rb") #The method tell() returns the current position of the file read/write pointer within the file.
276 #RecCounter=0
277 #Off2StartNxtRec=811248
278 #print 'OffsetStartHeader ',self.OffsetStartHeader,'RecCounter ', self.RecCounter, 'Off2StartNxtRec ' , self.Off2StartNxtRec
286 # print fp
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.
288 # startFp = open(fp,"rb") #The method tell() returns the current position of the file read/write pointer within the file.
289 # RecCounter=0
290 # Off2StartNxtRec=811248
291 # print 'OffsetStartHeader ',self.OffsetStartHeader,'RecCounter ', self.RecCounter, 'Off2StartNxtRec ' , self.Off2StartNxtRec
279 292 #OffRHeader= self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
280 293 #startFp.seek(OffRHeader, os.SEEK_SET)
281 print 'debe ser 48, RecCounter*811248', self.OffsetStartHeader,self.RecCounter,self.Off2StartNxtRec
282 print 'Posicion del bloque: ',OffRHeader
294 print 'debe ser 48, RecCounter*811248', self.OffsetStartHeader, self.RecCounter, self.Off2StartNxtRec
295 print 'Posicion del bloque: ', OffRHeader
283 296
284 header = numpy.fromfile(startFp,SRVI_STRUCTURE,1)
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] #
@@ -303,19 +316,16 self.northAngle = header['frame_cnt'][0] #
303 316 self.microsec = header['frame_cnt'][0] #
304 317 self.azisetvel = header['frame_cnt'][0] #
305 318 self.elvsetpos = header['frame_cnt'][0] #
306 self.RadarConst = header['frame_cnt'][0] #
319 self.RadarConst = header['frame_cnt'][0] #
307 320
308 321
309 self.ipp= 0.5*(SPEED_OF_LIGHT/self.PRFhz)
322 self.ipp = 0.5 * (SPEED_OF_LIGHT / self.PRFhz)
310 323
311 self.RHsize = 180+20*self.nChannels
312 self.Datasize= self.nProfiles*self.nChannels*self.nHeights*2*4
313 #print 'Datasize',self.Datasize
314 endFp = self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
324 self.RHsize = 180 + 20 * self.nChannels
325 self.Datasize = self.nProfiles * self.nChannels * self.nHeights * 2 * 4
326 # print 'Datasize',self.Datasize
327 endFp = self.OffsetStartHeader + self.RecCounter * self.Off2StartNxtRec
315 328
316 329 print '=============================================='
317 330
318 331 print '=============================================='
319
320
321 No newline at end of file
Show file before | Show file after | Hide comments
This diff has been collapsed as it changes many lines, (1270 lines changed) Show them Hide them
@@ -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
@@ -33,73 +33,71 from jroIO_base import JRODataReader
33 33
34 34
35 35 class Header(object):
36
36
37 37 def __init__(self):
38 38 raise NotImplementedError
39
40
39
41 40 def read(self):
42
41
43 42 raise NotImplementedError
44
43
45 44 def write(self):
46
45
47 46 raise NotImplementedError
48
47
49 48 def printInfo(self):
50
51 message = "#"*50 + "\n"
49
50 message = "#" * 50 + "\n"
52 51 message += self.__class__.__name__.upper() + "\n"
53 message += "#"*50 + "\n"
54
52 message += "#" * 50 + "\n"
53
55 54 keyList = self.__dict__.keys()
56 55 keyList.sort()
57
56
58 57 for key in keyList:
59 message += "%s = %s" %(key, self.__dict__[key]) + "\n"
60
58 message += "%s = %s" % (key, self.__dict__[key]) + "\n"
59
61 60 if "size" not in keyList:
62 61 attr = getattr(self, "size")
63
64 if attr:
65 message += "%s = %s" %("size", attr) + "\n"
66
67 #print message
68 62
63 if attr:
64 message += "%s = %s" % ("size", attr) + "\n"
69 65
66 # print message
70 67
71 68
69 FILE_STRUCTURE = numpy.dtype([ # HEADER 48bytes
70 ('FileMgcNumber', '<u4'), # 0x23020100
71 # No Of FDT data records in this file (0 or more)
72 ('nFDTdataRecors', '<u4'),
73 ('OffsetStartHeader', '<u4'),
74 ('RadarUnitId', '<u4'),
75 ('SiteName', numpy.str_, 32), # Null terminated
76 ])
72 77
73 FILE_STRUCTURE = numpy.dtype([ #HEADER 48bytes
74 ('FileMgcNumber','<u4'), #0x23020100
75 ('nFDTdataRecors','<u4'), #No Of FDT data records in this file (0 or more)
76 ('OffsetStartHeader','<u4'),
77 ('RadarUnitId','<u4'),
78 ('SiteName',numpy.str_,32), #Null terminated
79 ])
80 78
81 79 class FileHeaderBLTR(Header):
82
80
83 81 def __init__(self):
84
85 self.FileMgcNumber= 0 #0x23020100
86 self.nFDTdataRecors=0 #No Of FDT data records in this file (0 or more)
87 self.RadarUnitId= 0
88 self.OffsetStartHeader=0
89 self.SiteName= ""
82
83 self.FileMgcNumber = 0 # 0x23020100
84 # No Of FDT data records in this file (0 or more)
85 self.nFDTdataRecors = 0
86 self.RadarUnitId = 0
87 self.OffsetStartHeader = 0
88 self.SiteName = ""
90 89 self.size = 48
91
90
92 91 def FHread(self, fp):
93 #try:
94 startFp = open(fp,"rb")
95
96 header = numpy.fromfile(startFp, FILE_STRUCTURE,1)
97
92 # try:
93 startFp = open(fp, "rb")
94
95 header = numpy.fromfile(startFp, FILE_STRUCTURE, 1)
96
98 97 print ' '
99 98 print 'puntero file header', startFp.tell()
100 99 print ' '
101
102
100
103 101 ''' numpy.fromfile(file, dtype, count, sep='')
104 102 file : file or str
105 103 Open file object or filename.
@@ -119,32 +117,28 class FileHeaderBLTR(Header):
119 117
120 118 '''
121 119
120 self.FileMgcNumber = hex(header['FileMgcNumber'][0])
121 # No Of FDT data records in this file (0 or more)
122 self.nFDTdataRecors = int(header['nFDTdataRecors'][0])
123 self.RadarUnitId = int(header['RadarUnitId'][0])
124 self.OffsetStartHeader = int(header['OffsetStartHeader'][0])
125 self.SiteName = str(header['SiteName'][0])
126
127 # print 'Numero de bloques', self.nFDTdataRecors
122 128
123
124 self.FileMgcNumber= hex(header['FileMgcNumber'][0])
125 self.nFDTdataRecors=int(header['nFDTdataRecors'][0]) #No Of FDT data records in this file (0 or more)
126 self.RadarUnitId= int(header['RadarUnitId'][0])
127 self.OffsetStartHeader= int(header['OffsetStartHeader'][0])
128 self.SiteName= str(header['SiteName'][0])
129
130 #print 'Numero de bloques', self.nFDTdataRecors
131
132
133 if self.size <48:
129 if self.size < 48:
134 130 return 0
135
131
136 132 return 1
137
138
133
139 134 def write(self, fp):
140
135
141 136 headerTuple = (self.FileMgcNumber,
142 137 self.nFDTdataRecors,
143 138 self.RadarUnitId,
144 139 self.SiteName,
145 140 self.size)
146
147
141
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)
@@ -162,272 +156,308 class FileHeaderBLTR(Header):
162 156 first converting it to the closest Python type, and then using "format" % item.
163 157
164 158 '''
165
166 return 1
167 159
160 return 1
168 161
169 162
163 RECORD_STRUCTURE = numpy.dtype([ # RECORD HEADER 180+20N bytes
164 ('RecMgcNumber', '<u4'), # 0x23030001
165 ('RecCounter', '<u4'), # Record counter(0,1, ...)
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 ('TransmitFrec', '<u4'), # Transmit frequency (Hz)
185 ('ReceiveFrec', '<u4'), # Receive frequency
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'),
206 ('PRFhz', '<u4'), # PRF (Hz)
207 ('nCohInt', '<u4'), # Integrations
208 # Number of data points transformed
209 ('nProfiles', '<u4'),
210 # Number of receive beams stored in file (1 or N)
211 ('nChannels', '<u4'),
212 ('nIncohInt', '<u4'), # Number of spectral averages
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'),
249 ])
170 250
171 251
172 RECORD_STRUCTURE = numpy.dtype([ #RECORD HEADER 180+20N bytes
173 ('RecMgcNumber','<u4'), #0x23030001
174 ('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)
184 ('TransmitFrec','<u4'), #Transmit frequency (Hz)
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)
196 ('PRFhz','<u4'), #PRF (Hz)
197 ('nCohInt','<u4'), #Integrations
198 ('nProfiles','<u4'), #Number of data points transformed
199 ('nChannels','<u4'), #Number of receive beams stored in file (1 or N)
200 ('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
219 ])
252 class RecordHeaderBLTR(Header):
220 253
254 def __init__(self, RecMgcNumber=None, RecCounter=0, Off2StartNxtRec=811248,
255 nUtime=0, nMilisec=0, ExpTagName=None,
256 ExpComment=None, SiteLatDegrees=0, SiteLongDegrees=0,
257 RTCgpsStatus=0, TransmitFrec=0, ReceiveFrec=0,
258 FirstOsciFrec=0, Polarisation=0, ReceiverFiltSett=0,
259 nModesInUse=0, DualModeIndex=0, DualModeRange=0,
260 nDigChannels=0, SampResolution=0, nHeights=0,
261 StartRangeSamp=0, PRFhz=0, nCohInt=0,
262 nProfiles=0, nChannels=0, nIncohInt=0,
263 FFTwindowingInd=0, BeamAngleAzim=0, BeamAngleZen=0,
264 AntennaCoord0=0, AntennaCoord1=0, AntennaCoord2=0,
265 RecPhaseCalibr0=0, RecPhaseCalibr1=0, RecPhaseCalibr2=0,
266 RecAmpCalibr0=0, RecAmpCalibr1=0, RecAmpCalibr2=0,
267 AntennaAngl0=0, AntennaAngl1=0, AntennaAngl2=0,
268 ReceiverGaindB0=0, ReceiverGaindB1=0, ReceiverGaindB2=0, Off2StartData=0, OffsetStartHeader=0):
221 269
222 class RecordHeaderBLTR(Header):
223
224 def __init__(self, RecMgcNumber=None, RecCounter= 0, Off2StartNxtRec= 811248,
225 nUtime= 0, nMilisec= 0, ExpTagName= None,
226 ExpComment=None, SiteLatDegrees=0, SiteLongDegrees= 0,
227 RTCgpsStatus= 0, TransmitFrec= 0, ReceiveFrec= 0,
228 FirstOsciFrec= 0, Polarisation= 0, ReceiverFiltSett= 0,
229 nModesInUse= 0, DualModeIndex= 0, DualModeRange= 0,
230 nDigChannels= 0, SampResolution= 0, nHeights= 0,
231 StartRangeSamp= 0, PRFhz= 0, nCohInt= 0,
232 nProfiles= 0, nChannels= 0, nIncohInt= 0,
233 FFTwindowingInd= 0, BeamAngleAzim= 0, BeamAngleZen= 0,
234 AntennaCoord0= 0, AntennaCoord1= 0, AntennaCoord2= 0,
235 RecPhaseCalibr0= 0, RecPhaseCalibr1= 0, RecPhaseCalibr2= 0,
236 RecAmpCalibr0= 0, RecAmpCalibr1= 0, RecAmpCalibr2= 0,
237 AntennaAngl0=0, AntennaAngl1=0, AntennaAngl2=0,
238 ReceiverGaindB0= 0, ReceiverGaindB1= 0, ReceiverGaindB2= 0, Off2StartData=0, OffsetStartHeader=0):
239
240 self.RecMgcNumber = RecMgcNumber #0x23030001
241 self.RecCounter = RecCounter
242 self.Off2StartNxtRec = Off2StartNxtRec
270 self.RecMgcNumber = RecMgcNumber # 0x23030001
271 self.RecCounter = RecCounter
272 self.Off2StartNxtRec = Off2StartNxtRec
243 273 self.Off2StartData = Off2StartData
244 self.nUtime = nUtime
245 self.nMilisec = nMilisec
246 self.ExpTagName = ExpTagName
247 self.ExpComment = ExpComment
248 self.SiteLatDegrees = SiteLatDegrees
249 self.SiteLongDegrees = SiteLongDegrees
250 self.RTCgpsStatus = RTCgpsStatus
251 self.TransmitFrec = TransmitFrec
252 self.ReceiveFrec = ReceiveFrec
253 self.FirstOsciFrec = FirstOsciFrec
254 self.Polarisation = Polarisation
255 self.ReceiverFiltSett = ReceiverFiltSett
256 self.nModesInUse = nModesInUse
257 self.DualModeIndex = DualModeIndex
258 self.DualModeRange = DualModeRange
274 self.nUtime = nUtime
275 self.nMilisec = nMilisec
276 self.ExpTagName = ExpTagName
277 self.ExpComment = ExpComment
278 self.SiteLatDegrees = SiteLatDegrees
279 self.SiteLongDegrees = SiteLongDegrees
280 self.RTCgpsStatus = RTCgpsStatus
281 self.TransmitFrec = TransmitFrec
282 self.ReceiveFrec = ReceiveFrec
283 self.FirstOsciFrec = FirstOsciFrec
284 self.Polarisation = Polarisation
285 self.ReceiverFiltSett = ReceiverFiltSett
286 self.nModesInUse = nModesInUse
287 self.DualModeIndex = DualModeIndex
288 self.DualModeRange = DualModeRange
259 289 self.nDigChannels = nDigChannels
260 self.SampResolution = SampResolution
261 self.nHeights = nHeights
262 self.StartRangeSamp = StartRangeSamp
263 self.PRFhz = PRFhz
264 self.nCohInt = nCohInt
265 self.nProfiles = nProfiles
266 self.nChannels = nChannels
267 self.nIncohInt = nIncohInt
268 self.FFTwindowingInd = FFTwindowingInd
269 self.BeamAngleAzim = BeamAngleAzim
270 self.BeamAngleZen = BeamAngleZen
290 self.SampResolution = SampResolution
291 self.nHeights = nHeights
292 self.StartRangeSamp = StartRangeSamp
293 self.PRFhz = PRFhz
294 self.nCohInt = nCohInt
295 self.nProfiles = nProfiles
296 self.nChannels = nChannels
297 self.nIncohInt = nIncohInt
298 self.FFTwindowingInd = FFTwindowingInd
299 self.BeamAngleAzim = BeamAngleAzim
300 self.BeamAngleZen = BeamAngleZen
271 301 self.AntennaCoord0 = AntennaCoord0
272 302 self.AntennaAngl0 = AntennaAngl0
273 303 self.AntennaAngl1 = AntennaAngl1
274 304 self.AntennaAngl2 = AntennaAngl2
275 305 self.AntennaCoord1 = AntennaCoord1
276 self.AntennaCoord2 = AntennaCoord2
306 self.AntennaCoord2 = AntennaCoord2
277 307 self.RecPhaseCalibr0 = RecPhaseCalibr0
278 308 self.RecPhaseCalibr1 = RecPhaseCalibr1
279 self.RecPhaseCalibr2 = RecPhaseCalibr2
309 self.RecPhaseCalibr2 = RecPhaseCalibr2
280 310 self.RecAmpCalibr0 = RecAmpCalibr0
281 311 self.RecAmpCalibr1 = RecAmpCalibr1
282 self.RecAmpCalibr2 = RecAmpCalibr2
312 self.RecAmpCalibr2 = RecAmpCalibr2
283 313 self.ReceiverGaindB0 = ReceiverGaindB0
284 314 self.ReceiverGaindB1 = ReceiverGaindB1
285 self.ReceiverGaindB2 = ReceiverGaindB2
286 self.OffsetStartHeader = 48
287
288
289
315 self.ReceiverGaindB2 = ReceiverGaindB2
316 self.OffsetStartHeader = 48
317
290 318 def RHread(self, fp):
291 #print fp
292 #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.
294 #RecCounter=0
295 #Off2StartNxtRec=811248
296 OffRHeader= self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
319 # print fp
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.
321 # The method tell() returns the current position of the file read/write pointer within the file.
322 startFp = open(fp, "rb")
323 # RecCounter=0
324 # Off2StartNxtRec=811248
325 OffRHeader = self.OffsetStartHeader + self.RecCounter * self.Off2StartNxtRec
297 326 print ' '
298 327 print 'puntero Record Header', startFp.tell()
299 328 print ' '
300
301
329
302 330 startFp.seek(OffRHeader, os.SEEK_SET)
303
331
304 332 print ' '
305 333 print 'puntero Record Header con seek', startFp.tell()
306 334 print ' '
307
308 #print 'Posicion del bloque: ',OffRHeader
309
310 header = numpy.fromfile(startFp,RECORD_STRUCTURE,1)
311
335
336 # print 'Posicion del bloque: ',OffRHeader
337
338 header = numpy.fromfile(startFp, RECORD_STRUCTURE, 1)
339
312 340 print ' '
313 341 print 'puntero Record Header con seek', startFp.tell()
314 342 print ' '
315
343
316 344 print ' '
317 345 #
318 #print 'puntero Record Header despues de seek', header.tell()
346 # print 'puntero Record Header despues de seek', header.tell()
319 347 print ' '
320
321 self.RecMgcNumber = hex(header['RecMgcNumber'][0]) #0x23030001
322 self.RecCounter = int(header['RecCounter'][0])
323 self.Off2StartNxtRec = int(header['Off2StartNxtRec'][0])
324 self.Off2StartData = int(header['Off2StartData'][0])
325 self.nUtime = header['nUtime'][0]
326 self.nMilisec = header['nMilisec'][0]
327 self.ExpTagName = str(header['ExpTagName'][0])
328 self.ExpComment = str(header['ExpComment'][0])
329 self.SiteLatDegrees = header['SiteLatDegrees'][0]
330 self.SiteLongDegrees = header['SiteLongDegrees'][0]
331 self.RTCgpsStatus = header['RTCgpsStatus'][0]
332 self.TransmitFrec = header['TransmitFrec'][0]
333 self.ReceiveFrec = header['ReceiveFrec'][0]
334 self.FirstOsciFrec = header['FirstOsciFrec'][0]
335 self.Polarisation = header['Polarisation'][0]
336 self.ReceiverFiltSett = header['ReceiverFiltSett'][0]
337 self.nModesInUse = header['nModesInUse'][0]
338 self.DualModeIndex = header['DualModeIndex'][0]
339 self.DualModeRange = header['DualModeRange'][0]
348
349 self.RecMgcNumber = hex(header['RecMgcNumber'][0]) # 0x23030001
350 self.RecCounter = int(header['RecCounter'][0])
351 self.Off2StartNxtRec = int(header['Off2StartNxtRec'][0])
352 self.Off2StartData = int(header['Off2StartData'][0])
353 self.nUtime = header['nUtime'][0]
354 self.nMilisec = header['nMilisec'][0]
355 self.ExpTagName = str(header['ExpTagName'][0])
356 self.ExpComment = str(header['ExpComment'][0])
357 self.SiteLatDegrees = header['SiteLatDegrees'][0]
358 self.SiteLongDegrees = header['SiteLongDegrees'][0]
359 self.RTCgpsStatus = header['RTCgpsStatus'][0]
360 self.TransmitFrec = header['TransmitFrec'][0]
361 self.ReceiveFrec = header['ReceiveFrec'][0]
362 self.FirstOsciFrec = header['FirstOsciFrec'][0]
363 self.Polarisation = header['Polarisation'][0]
364 self.ReceiverFiltSett = header['ReceiverFiltSett'][0]
365 self.nModesInUse = header['nModesInUse'][0]
366 self.DualModeIndex = header['DualModeIndex'][0]
367 self.DualModeRange = header['DualModeRange'][0]
340 368 self.nDigChannels = header['nDigChannels'][0]
341 self.SampResolution = header['SampResolution'][0]
342 self.nHeights = header['nHeights'][0]
343 self.StartRangeSamp = header['StartRangeSamp'][0]
344 self.PRFhz = header['PRFhz'][0]
345 self.nCohInt = header['nCohInt'][0]
346 self.nProfiles = header['nProfiles'][0]
347 self.nChannels = header['nChannels'][0]
348 self.nIncohInt = header['nIncohInt'][0]
349 self.FFTwindowingInd = header['FFTwindowingInd'][0]
350 self.BeamAngleAzim = header['BeamAngleAzim'][0]
351 self.BeamAngleZen = header['BeamAngleZen'][0]
369 self.SampResolution = header['SampResolution'][0]
370 self.nHeights = header['nHeights'][0]
371 self.StartRangeSamp = header['StartRangeSamp'][0]
372 self.PRFhz = header['PRFhz'][0]
373 self.nCohInt = header['nCohInt'][0]
374 self.nProfiles = header['nProfiles'][0]
375 self.nChannels = header['nChannels'][0]
376 self.nIncohInt = header['nIncohInt'][0]
377 self.FFTwindowingInd = header['FFTwindowingInd'][0]
378 self.BeamAngleAzim = header['BeamAngleAzim'][0]
379 self.BeamAngleZen = header['BeamAngleZen'][0]
352 380 self.AntennaCoord0 = header['AntennaCoord0'][0]
353 381 self.AntennaAngl0 = header['AntennaAngl0'][0]
354 382 self.AntennaCoord1 = header['AntennaCoord1'][0]
355 self.AntennaAngl1 = header['AntennaAngl1'][0]
383 self.AntennaAngl1 = header['AntennaAngl1'][0]
356 384 self.AntennaCoord2 = header['AntennaCoord2'][0]
357 self.AntennaAngl2 = header['AntennaAngl2'][0]
385 self.AntennaAngl2 = header['AntennaAngl2'][0]
358 386 self.RecPhaseCalibr0 = header['RecPhaseCalibr0'][0]
359 387 self.RecPhaseCalibr1 = header['RecPhaseCalibr1'][0]
360 self.RecPhaseCalibr2 = header['RecPhaseCalibr2'][0]
388 self.RecPhaseCalibr2 = header['RecPhaseCalibr2'][0]
361 389 self.RecAmpCalibr0 = header['RecAmpCalibr0'][0]
362 390 self.RecAmpCalibr1 = header['RecAmpCalibr1'][0]
363 self.RecAmpCalibr2 = header['RecAmpCalibr2'][0]
391 self.RecAmpCalibr2 = header['RecAmpCalibr2'][0]
364 392 self.ReceiverGaindB0 = header['ReceiverGaindB0'][0]
365 393 self.ReceiverGaindB1 = header['ReceiverGaindB1'][0]
366 394 self.ReceiverGaindB2 = header['ReceiverGaindB2'][0]
367 395
368 self.ipp= 0.5*(SPEED_OF_LIGHT/self.PRFhz)
369
370 self.RHsize = 180+20*self.nChannels
371 self.Datasize= self.nProfiles*self.nChannels*self.nHeights*2*4
372 #print 'Datasize',self.Datasize
373 endFp = self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
374
396 self.ipp = 0.5 * (SPEED_OF_LIGHT / self.PRFhz)
397
398 self.RHsize = 180 + 20 * self.nChannels
399 self.Datasize = self.nProfiles * self.nChannels * self.nHeights * 2 * 4
400 # print 'Datasize',self.Datasize
401 endFp = self.OffsetStartHeader + self.RecCounter * self.Off2StartNxtRec
402
375 403 print '=============================================='
376 print 'RecMgcNumber ',self.RecMgcNumber
377 print 'RecCounter ',self.RecCounter
378 print 'Off2StartNxtRec ',self.Off2StartNxtRec
379 print 'Off2StartData ',self.Off2StartData
380 print 'Range Resolution ',self.SampResolution
381 print 'First Height ',self.StartRangeSamp
382 print 'PRF (Hz) ',self.PRFhz
383 print 'Heights (K) ',self.nHeights
384 print 'Channels (N) ',self.nChannels
385 print 'Profiles (J) ',self.nProfiles
386 print 'iCoh ',self.nCohInt
387 print 'iInCoh ',self.nIncohInt
388 print 'BeamAngleAzim ',self.BeamAngleAzim
389 print 'BeamAngleZen ',self.BeamAngleZen
390
391 #print 'ModoEnUso ',self.DualModeIndex
392 #print 'UtcTime ',self.nUtime
393 #print 'MiliSec ',self.nMilisec
394 #print 'Exp TagName ',self.ExpTagName
395 #print 'Exp Comment ',self.ExpComment
396 #print 'FFT Window Index ',self.FFTwindowingInd
397 #print 'N Dig. Channels ',self.nDigChannels
398 print 'Size de bloque ',self.RHsize
399 print 'DataSize ',self.Datasize
400 print 'BeamAngleAzim ',self.BeamAngleAzim
401 #print 'AntennaCoord0 ',self.AntennaCoord0
402 #print 'AntennaAngl0 ',self.AntennaAngl0
403 #print 'AntennaCoord1 ',self.AntennaCoord1
404 #print 'AntennaAngl1 ',self.AntennaAngl1
405 #print 'AntennaCoord2 ',self.AntennaCoord2
406 #print 'AntennaAngl2 ',self.AntennaAngl2
407 print 'RecPhaseCalibr0 ',self.RecPhaseCalibr0
408 print 'RecPhaseCalibr1 ',self.RecPhaseCalibr1
409 print 'RecPhaseCalibr2 ',self.RecPhaseCalibr2
410 print 'RecAmpCalibr0 ',self.RecAmpCalibr0
411 print 'RecAmpCalibr1 ',self.RecAmpCalibr1
412 print 'RecAmpCalibr2 ',self.RecAmpCalibr2
413 print 'ReceiverGaindB0 ',self.ReceiverGaindB0
414 print 'ReceiverGaindB1 ',self.ReceiverGaindB1
415 print 'ReceiverGaindB2 ',self.ReceiverGaindB2
404 print 'RecMgcNumber ', self.RecMgcNumber
405 print 'RecCounter ', self.RecCounter
406 print 'Off2StartNxtRec ', self.Off2StartNxtRec
407 print 'Off2StartData ', self.Off2StartData
408 print 'Range Resolution ', self.SampResolution
409 print 'First Height ', self.StartRangeSamp
410 print 'PRF (Hz) ', self.PRFhz
411 print 'Heights (K) ', self.nHeights
412 print 'Channels (N) ', self.nChannels
413 print 'Profiles (J) ', self.nProfiles
414 print 'iCoh ', self.nCohInt
415 print 'iInCoh ', self.nIncohInt
416 print 'BeamAngleAzim ', self.BeamAngleAzim
417 print 'BeamAngleZen ', self.BeamAngleZen
418
419 # print 'ModoEnUso ',self.DualModeIndex
420 # print 'UtcTime ',self.nUtime
421 # print 'MiliSec ',self.nMilisec
422 # print 'Exp TagName ',self.ExpTagName
423 # print 'Exp Comment ',self.ExpComment
424 # print 'FFT Window Index ',self.FFTwindowingInd
425 # print 'N Dig. Channels ',self.nDigChannels
426 print 'Size de bloque ', self.RHsize
427 print 'DataSize ', self.Datasize
428 print 'BeamAngleAzim ', self.BeamAngleAzim
429 # print 'AntennaCoord0 ',self.AntennaCoord0
430 # print 'AntennaAngl0 ',self.AntennaAngl0
431 # print 'AntennaCoord1 ',self.AntennaCoord1
432 # print 'AntennaAngl1 ',self.AntennaAngl1
433 # print 'AntennaCoord2 ',self.AntennaCoord2
434 # print 'AntennaAngl2 ',self.AntennaAngl2
435 print 'RecPhaseCalibr0 ', self.RecPhaseCalibr0
436 print 'RecPhaseCalibr1 ', self.RecPhaseCalibr1
437 print 'RecPhaseCalibr2 ', self.RecPhaseCalibr2
438 print 'RecAmpCalibr0 ', self.RecAmpCalibr0
439 print 'RecAmpCalibr1 ', self.RecAmpCalibr1
440 print 'RecAmpCalibr2 ', self.RecAmpCalibr2
441 print 'ReceiverGaindB0 ', self.ReceiverGaindB0
442 print 'ReceiverGaindB1 ', self.ReceiverGaindB1
443 print 'ReceiverGaindB2 ', self.ReceiverGaindB2
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
427 457
428
458
429 459 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODataReader):
430
460
431 461 path = None
432 462 startDate = None
433 463 endDate = None
@@ -435,28 +465,25 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
435 465 endTime = None
436 466 walk = None
437 467 isConfig = False
438
439
440 fileList= None
441
442 #metadata
443 TimeZone= None
444 Interval= None
445 heightList= None
446
447 #data
448 data= None
449 utctime= None
450
451
452
468
469 fileList = None
470
471 # metadata
472 TimeZone = None
473 Interval = None
474 heightList = None
475
476 # data
477 data = None
478 utctime = None
479
453 480 def __init__(self, **kwargs):
454
455 #Eliminar de la base la herencia
481
482 # Eliminar de la base la herencia
456 483 ProcessingUnit.__init__(self, **kwargs)
457
484
458 485 #self.isConfig = False
459
486
460 487 #self.pts2read_SelfSpectra = 0
461 488 #self.pts2read_CrossSpectra = 0
462 489 #self.pts2read_DCchannels = 0
@@ -464,60 +491,59 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
464 491 self.utc = None
465 492 self.ext = ".fdt"
466 493 self.optchar = "P"
467 self.fpFile=None
494 self.fpFile = None
468 495 self.fp = None
469 self.BlockCounter=0
496 self.BlockCounter = 0
470 497 self.dtype = None
471 498 self.fileSizeByHeader = None
472 499 self.filenameList = []
473 500 self.fileSelector = 0
474 self.Off2StartNxtRec=0
475 self.RecCounter=0
501 self.Off2StartNxtRec = 0
502 self.RecCounter = 0
476 503 self.flagNoMoreFiles = 0
477 self.data_spc=None
478 self.data_cspc=None
479 self.data_output=None
504 self.data_spc = None
505 self.data_cspc = None
506 self.data_output = None
480 507 self.path = None
481 self.OffsetStartHeader=0
482 self.Off2StartData=0
508 self.OffsetStartHeader = 0
509 self.Off2StartData = 0
483 510 self.ipp = 0
484 self.nFDTdataRecors=0
511 self.nFDTdataRecors = 0
485 512 self.blocksize = 0
486 513 self.dataOut = Spectra()
487 self.profileIndex = 1 #Always
488 self.dataOut.flagNoData=False
514 self.profileIndex = 1 # Always
515 self.dataOut.flagNoData = False
489 516 self.dataOut.nRdPairs = 0
490 517 self.dataOut.pairsList = []
491 self.dataOut.data_spc=None
492 self.dataOut.noise=[]
493 self.dataOut.velocityX=[]
494 self.dataOut.velocityY=[]
495 self.dataOut.velocityV=[]
496
497
518 self.dataOut.data_spc = None
519 self.dataOut.noise = []
520 self.dataOut.velocityX = []
521 self.dataOut.velocityY = []
522 self.dataOut.velocityV = []
498 523
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.
502 527 It also creates an organized list with the names of the files to read.
503 528 '''
504 #self.__checkPath()
505
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
508 nFiles=0 #File Counter
509 FileList=[] #A list is created that will contain the .fdt files
510 for IndexFile in ListaData :
511 if '.fdt' in IndexFile:
529 # self.__checkPath()
530
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)
535 nFiles = 0 # File Counter
536 FileList = [] # A list is created that will contain the .fdt files
537 for IndexFile in ListaData:
538 if '.fdt' in IndexFile:
512 539 FileList.append(IndexFile)
513 nFiles+=1
514
515 #print 'Files2Read'
516 #print 'Existen '+str(nFiles)+' archivos .fdt'
517
518 self.filenameList=FileList #List of files from least to largest by names
519
520
540 nFiles += 1
541
542 # print 'Files2Read'
543 # print 'Existen '+str(nFiles)+' archivos .fdt'
544
545 self.filenameList = FileList # List of files from least to largest by names
546
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.
@@ -527,341 +553,350 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
527 553 if not self.isConfig:
528 554 self.setup(**kwargs)
529 555 self.isConfig = True
530
556
531 557 self.getData()
532 #print 'running'
533
534
558 # print 'running'
559
535 560 def setup(self, path=None,
536 startDate=None,
537 endDate=None,
538 startTime=None,
539 endTime=None,
540 walk=True,
541 timezone='utc',
542 code = None,
543 online=False,
544 ReadMode=None,
545 **kwargs):
546
561 startDate=None,
562 endDate=None,
563 startTime=None,
564 endTime=None,
565 walk=True,
566 timezone='utc',
567 code=None,
568 online=False,
569 ReadMode=None,
570 **kwargs):
571
547 572 self.isConfig = True
548
549 self.path=path
550 self.startDate=startDate
551 self.endDate=endDate
552 self.startTime=startTime
553 self.endTime=endTime
554 self.walk=walk
555 self.ReadMode=int(ReadMode)
556
573
574 self.path = path
575 self.startDate = startDate
576 self.endDate = endDate
577 self.startTime = startTime
578 self.endTime = endTime
579 self.walk = walk
580 self.ReadMode = int(ReadMode)
581
557 582 pass
558
559
583
560 584 def getData(self):
561 585 '''
562 586 Before starting this function, you should check that there is still an unread file,
563 587 If there are still blocks to read or if the data block is empty.
564
588
565 589 You should call the file "read".
566
590
567 591 '''
568
592
569 593 if self.flagNoMoreFiles:
570 594 self.dataOut.flagNoData = True
571 595 print 'NoData se vuelve true'
572 596 return 0
573 597
574 self.fp=self.path
598 self.fp = self.path
575 599 self.Files2Read(self.fp)
576 600 self.readFile(self.fp)
577 601 self.dataOut.data_spc = self.data_spc
578 self.dataOut.data_cspc =self.data_cspc
579 self.dataOut.data_output=self.data_output
580
602 self.dataOut.data_cspc = self.data_cspc
603 self.dataOut.data_output = self.data_output
604
581 605 print 'self.dataOut.data_output', shape(self.dataOut.data_output)
582
583 #self.removeDC()
584 return self.dataOut.data_spc
585
586
587 def readFile(self,fp):
606
607 # self.removeDC()
608 return self.dataOut.data_spc
609
610 def readFile(self, fp):
588 611 '''
589 612 You must indicate if you are reading in Online or Offline mode and load the
590 613 The parameters for this file reading mode.
591
614
592 615 Then you must do 2 actions:
593
616
594 617 1. Get the BLTR FileHeader.
595 618 2. Start reading the first block.
596 619 '''
597
598 #The address of the folder is generated the name of the .fdt file that will be read
599 print "File: ",self.fileSelector+1
600
620
621 # The address of the folder is generated the name of the .fdt file that will be read
622 print "File: ", self.fileSelector + 1
623
601 624 if self.fileSelector < len(self.filenameList):
602
603 self.fpFile=str(fp)+'/'+str(self.filenameList[self.fileSelector])
604 #print self.fpFile
625
626 self.fpFile = str(fp) + '/' + \
627 str(self.filenameList[self.fileSelector])
628 # print self.fpFile
605 629 fheader = FileHeaderBLTR()
606 fheader.FHread(self.fpFile) #Bltr FileHeader Reading
607 self.nFDTdataRecors=fheader.nFDTdataRecors
608
609 self.readBlock() #Block reading
630 fheader.FHread(self.fpFile) # Bltr FileHeader Reading
631 self.nFDTdataRecors = fheader.nFDTdataRecors
632
633 self.readBlock() # Block reading
610 634 else:
611 635 print 'readFile FlagNoData becomes true'
612 self.flagNoMoreFiles=True
636 self.flagNoMoreFiles = True
613 637 self.dataOut.flagNoData = True
614 return 0
615
638 return 0
639
616 640 def getVelRange(self, extrapoints=0):
617 Lambda= SPEED_OF_LIGHT/50000000
618 PRF = self.dataOut.PRF#1./(self.dataOut.ippSeconds * self.dataOut.nCohInt)
619 Vmax=-Lambda/(4.*(1./PRF)*self.dataOut.nCohInt*2.)
620 deltafreq = PRF / (self.nProfiles)
621 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.)
624 return velrange
625
626 def readBlock(self):
641 Lambda = SPEED_OF_LIGHT / 50000000
642 # 1./(self.dataOut.ippSeconds * self.dataOut.nCohInt)
643 PRF = self.dataOut.PRF
644 Vmax = -Lambda / (4. * (1. / PRF) * self.dataOut.nCohInt * 2.)
645 deltafreq = PRF / (self.nProfiles)
646 deltavel = (Vmax * 2) / (self.nProfiles)
647 freqrange = deltafreq * \
648 (numpy.arange(self.nProfiles) - self.nProfiles / 2.) - deltafreq / 2
649 velrange = deltavel * \
650 (numpy.arange(self.nProfiles) - self.nProfiles / 2.)
651 return velrange
652
653 def readBlock(self):
627 654 '''
628 655 It should be checked if the block has data, if it is not passed to the next file.
629
656
630 657 Then the following is done:
631
658
632 659 1. Read the RecordHeader
633 660 2. Fill the buffer with the current block number.
634
661
635 662 '''
636
637 if self.BlockCounter < self.nFDTdataRecors-2:
663
664 if self.BlockCounter < self.nFDTdataRecors - 2:
638 665 print self.nFDTdataRecors, 'CONDICION!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
639 if self.ReadMode==1:
640 rheader = RecordHeaderBLTR(RecCounter=self.BlockCounter+1)
641 elif self.ReadMode==0:
666 if self.ReadMode == 1:
667 rheader = RecordHeaderBLTR(RecCounter=self.BlockCounter + 1)
668 elif self.ReadMode == 0:
642 669 rheader = RecordHeaderBLTR(RecCounter=self.BlockCounter)
643
644 rheader.RHread(self.fpFile) #Bltr FileHeader Reading
645
646 self.OffsetStartHeader=rheader.OffsetStartHeader
647 self.RecCounter=rheader.RecCounter
648 self.Off2StartNxtRec=rheader.Off2StartNxtRec
649 self.Off2StartData=rheader.Off2StartData
650 self.nProfiles=rheader.nProfiles
651 self.nChannels=rheader.nChannels
652 self.nHeights=rheader.nHeights
653 self.frequency=rheader.TransmitFrec
654 self.DualModeIndex=rheader.DualModeIndex
655
656 self.pairsList =[(0,1),(0,2),(1,2)]
670
671 rheader.RHread(self.fpFile) # Bltr FileHeader Reading
672
673 self.OffsetStartHeader = rheader.OffsetStartHeader
674 self.RecCounter = rheader.RecCounter
675 self.Off2StartNxtRec = rheader.Off2StartNxtRec
676 self.Off2StartData = rheader.Off2StartData
677 self.nProfiles = rheader.nProfiles
678 self.nChannels = rheader.nChannels
679 self.nHeights = rheader.nHeights
680 self.frequency = rheader.TransmitFrec
681 self.DualModeIndex = rheader.DualModeIndex
682
683 self.pairsList = [(0, 1), (0, 2), (1, 2)]
657 684 self.dataOut.pairsList = self.pairsList
658
659 self.nRdPairs=len(self.dataOut.pairsList)
685
686 self.nRdPairs = len(self.dataOut.pairsList)
660 687 self.dataOut.nRdPairs = self.nRdPairs
661
662 self.__firstHeigth=rheader.StartRangeSamp
663 self.__deltaHeigth=rheader.SampResolution
664 self.dataOut.heightList= self.__firstHeigth + numpy.array(range(self.nHeights))*self.__deltaHeigth
688
689 self.__firstHeigth = rheader.StartRangeSamp
690 self.__deltaHeigth = rheader.SampResolution
691 self.dataOut.heightList = self.__firstHeigth + \
692 numpy.array(range(self.nHeights)) * self.__deltaHeigth
665 693 self.dataOut.channelList = range(self.nChannels)
666 self.dataOut.nProfiles=rheader.nProfiles
667 self.dataOut.nIncohInt=rheader.nIncohInt
668 self.dataOut.nCohInt=rheader.nCohInt
669 self.dataOut.ippSeconds= 1/float(rheader.PRFhz)
670 self.dataOut.PRF=rheader.PRFhz
671 self.dataOut.nFFTPoints=rheader.nProfiles
672 self.dataOut.utctime=rheader.nUtime
673 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
676
677 self.data_output=numpy.ones([3,rheader.nHeights])*numpy.NaN
694 self.dataOut.nProfiles = rheader.nProfiles
695 self.dataOut.nIncohInt = rheader.nIncohInt
696 self.dataOut.nCohInt = rheader.nCohInt
697 self.dataOut.ippSeconds = 1 / float(rheader.PRFhz)
698 self.dataOut.PRF = rheader.PRFhz
699 self.dataOut.nFFTPoints = rheader.nProfiles
700 self.dataOut.utctime = rheader.nUtime
701 self.dataOut.timeZone = 0
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
706
707 self.data_output = numpy.ones([3, rheader.nHeights]) * numpy.NaN
678 708 print 'self.data_output', shape(self.data_output)
679 self.dataOut.velocityX=[]
680 self.dataOut.velocityY=[]
681 self.dataOut.velocityV=[]
682
709 self.dataOut.velocityX = []
710 self.dataOut.velocityY = []
711 self.dataOut.velocityV = []
712
683 713 '''Block Reading, the Block Data is received and Reshape is used to give it
684 714 shape.
685 '''
686
687 #Procedure to take the pointer to where the date block starts
688 startDATA = open(self.fpFile,"rb")
689 OffDATA= self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec+self.Off2StartData
715 '''
716
717 # Procedure to take the pointer to where the date block starts
718 startDATA = open(self.fpFile, "rb")
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):
693 return numpy.convolve(x, numpy.ones((N,))/N)[(N-1):]
694
695 def gaus(xSamples,a,x0,sigma):
696 return a*exp(-(xSamples-x0)**2/(2*sigma**2))
697
698 def Find(x,value):
724 return numpy.convolve(x, numpy.ones((N,)) / N)[(N - 1):]
725
726 def gaus(xSamples, a, x0, sigma):
727 return a * exp(-(xSamples - x0)**2 / (2 * sigma**2))
728
729 def Find(x, value):
699 730 for index in range(len(x)):
700 if x[index]==value:
701 return index
702
731 if x[index] == value:
732 return index
733
703 734 def pol2cart(rho, phi):
704 735 x = rho * numpy.cos(phi)
705 736 y = rho * numpy.sin(phi)
706 737 return(x, y)
707
708
709
710
711 if self.DualModeIndex==self.ReadMode:
712
713 self.data_fft = numpy.fromfile( startDATA, [('complex','<c8')],self.nProfiles*self.nChannels*self.nHeights )
714
715 self.data_fft=self.data_fft.astype(numpy.dtype('complex'))
716
717 self.data_block=numpy.reshape(self.data_fft,(self.nHeights, self.nChannels, self.nProfiles ))
718
719 self.data_block = numpy.transpose(self.data_block, (1,2,0))
720
738
739 if self.DualModeIndex == self.ReadMode:
740
741 self.data_fft = numpy.fromfile(
742 startDATA, [('complex', '<c8')], self.nProfiles * self.nChannels * self.nHeights)
743
744 self.data_fft = self.data_fft.astype(numpy.dtype('complex'))
745
746 self.data_block = numpy.reshape(
747 self.data_fft, (self.nHeights, self.nChannels, self.nProfiles))
748
749 self.data_block = numpy.transpose(self.data_block, (1, 2, 0))
750
721 751 copy = self.data_block.copy()
722 spc = copy * numpy.conjugate(copy)
723
724 self.data_spc = numpy.absolute(spc) # valor absoluto o magnitud
725
752 spc = copy * numpy.conjugate(copy)
753
754 self.data_spc = numpy.absolute(
755 spc) # valor absoluto o magnitud
756
726 757 factor = self.dataOut.normFactor
727
728
729 z = self.data_spc.copy()#/factor
758
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)
732 762 print ' '
733 763 print 'Z: '
734 print shape(z)
764 print shape(z)
735 765 print ' '
736 766 print ' '
737
738 self.dataOut.data_spc=self.data_spc
739
740 self.noise = self.dataOut.getNoise(ymin_index=80, ymax_index=132)#/factor
767
768 self.dataOut.data_spc = self.data_spc
769
770 self.noise = self.dataOut.getNoise(
771 ymin_index=80, ymax_index=132) # /factor
741 772 #noisedB = 10*numpy.log10(self.noise)
742
743
744 ySamples=numpy.ones([3,self.nProfiles])
745 phase=numpy.ones([3,self.nProfiles])
746 CSPCSamples=numpy.ones([3,self.nProfiles],dtype=numpy.complex_)
747 coherence=numpy.ones([3,self.nProfiles])
748 PhaseSlope=numpy.ones(3)
749 PhaseInter=numpy.ones(3)
750
773
774 ySamples = numpy.ones([3, self.nProfiles])
775 phase = numpy.ones([3, self.nProfiles])
776 CSPCSamples = numpy.ones(
777 [3, self.nProfiles], dtype=numpy.complex_)
778 coherence = numpy.ones([3, self.nProfiles])
779 PhaseSlope = numpy.ones(3)
780 PhaseInter = numpy.ones(3)
781
751 782 '''****** Getting CrossSpectra ******'''
752 cspc=self.data_block.copy()
753 self.data_cspc=self.data_block.copy()
754
755 xFrec=self.getVelRange(1)
756 VelRange=self.getVelRange(1)
757 self.dataOut.VelRange=VelRange
758 #print ' '
759 #print ' '
760 #print 'xFrec',xFrec
761 #print ' '
762 #print ' '
763 #Height=35
764 for i in range(self.nRdPairs):
765
783 cspc = self.data_block.copy()
784 self.data_cspc = self.data_block.copy()
785
786 xFrec = self.getVelRange(1)
787 VelRange = self.getVelRange(1)
788 self.dataOut.VelRange = VelRange
789 # print ' '
790 # print ' '
791 # print 'xFrec',xFrec
792 # print ' '
793 # print ' '
794 # Height=35
795 for i in range(self.nRdPairs):
796
766 797 chan_index0 = self.dataOut.pairsList[i][0]
767 798 chan_index1 = self.dataOut.pairsList[i][1]
768
769 self.data_cspc[i,:,:]=cspc[chan_index0,:,:] * numpy.conjugate(cspc[chan_index1,:,:])
770
771
799
800 self.data_cspc[i, :, :] = cspc[chan_index0, :,
801 :] * numpy.conjugate(cspc[chan_index1, :, :])
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)
776
777 E01=AntennaX0-AntennaX1
778 N01=AntennaY0-AntennaY1
779
780 E02=AntennaX0-AntennaX2
781 N02=AntennaY0-AntennaY2
782
783 E12=AntennaX1-AntennaX2
784 N12=AntennaY1-AntennaY2
785
786 self.ChanDist= numpy.array([[E01, N01],[E02,N02],[E12,N12]])
787
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)
810
811 E01 = AntennaX0 - AntennaX1
812 N01 = AntennaY0 - AntennaY1
813
814 E02 = AntennaX0 - AntennaX2
815 N02 = AntennaY0 - AntennaY2
816
817 E12 = AntennaX1 - AntennaX2
818 N12 = AntennaY1 - AntennaY2
819
820 self.ChanDist = numpy.array(
821 [[E01, N01], [E02, N02], [E12, N12]])
822
788 823 self.dataOut.ChanDist = self.ChanDist
789
790
824
825
791 826 # for Height in range(self.nHeights):
792 #
827 #
793 828 # for i in range(self.nRdPairs):
794 #
829 #
795 830 # '''****** Line of Data SPC ******'''
796 831 # zline=z[i,:,Height]
797 #
832 #
798 833 # '''****** DC is removed ******'''
799 834 # DC=Find(zline,numpy.amax(zline))
800 835 # zline[DC]=(zline[DC-1]+zline[DC+1])/2
801 #
802 #
836 #
837 #
803 838 # '''****** SPC is normalized ******'''
804 839 # FactNorm= zline.copy() / numpy.sum(zline.copy())
805 840 # FactNorm= FactNorm/numpy.sum(FactNorm)
806 #
841 #
807 842 # SmoothSPC=moving_average(FactNorm,N=3)
808 #
843 #
809 844 # xSamples = ar(range(len(SmoothSPC)))
810 845 # ySamples[i] = SmoothSPC-self.noise[i]
811 #
846 #
812 847 # for i in range(self.nRdPairs):
813 #
848 #
814 849 # '''****** Line of Data CSPC ******'''
815 850 # cspcLine=self.data_cspc[i,:,Height].copy()
816 #
817 #
818 #
851 #
852 #
853 #
819 854 # '''****** CSPC is normalized ******'''
820 855 # chan_index0 = self.dataOut.pairsList[i][0]
821 # chan_index1 = self.dataOut.pairsList[i][1]
856 # chan_index1 = self.dataOut.pairsList[i][1]
822 857 # CSPCFactor= numpy.sum(ySamples[chan_index0]) * numpy.sum(ySamples[chan_index1])
823 #
824 #
858 #
859 #
825 860 # CSPCNorm= cspcLine.copy() / numpy.sqrt(CSPCFactor)
826 #
827 #
861 #
862 #
828 863 # CSPCSamples[i] = CSPCNorm-self.noise[i]
829 864 # coherence[i] = numpy.abs(CSPCSamples[i]) / numpy.sqrt(CSPCFactor)
830 #
865 #
831 866 # '''****** DC is removed ******'''
832 867 # DC=Find(coherence[i],numpy.amax(coherence[i]))
833 868 # coherence[i][DC]=(coherence[i][DC-1]+coherence[i][DC+1])/2
834 869 # coherence[i]= moving_average(coherence[i],N=2)
835 #
870 #
836 871 # phase[i] = moving_average( numpy.arctan2(CSPCSamples[i].imag, CSPCSamples[i].real),N=1)#*180/numpy.pi
837 #
838 #
872 #
873 #
839 874 # '''****** Getting fij width ******'''
840 #
875 #
841 876 # yMean=[]
842 # yMean2=[]
843 #
877 # yMean2=[]
878 #
844 879 # for j in range(len(ySamples[1])):
845 880 # yMean=numpy.append(yMean,numpy.average([ySamples[0,j],ySamples[1,j],ySamples[2,j]]))
846 #
881 #
847 882 # '''******* Getting fitting Gaussian ******'''
848 883 # meanGauss=sum(xSamples*yMean) / len(xSamples)
849 884 # sigma=sum(yMean*(xSamples-meanGauss)**2) / len(xSamples)
850 885 # #print 'Height',Height,'SNR', meanGauss/sigma**2
851 #
886 #
852 887 # if (abs(meanGauss/sigma**2) > 0.0001) :
853 #
854 # try:
888 #
889 # try:
855 890 # popt,pcov = curve_fit(gaus,xSamples,yMean,p0=[1,meanGauss,sigma])
856 #
891 #
857 892 # if numpy.amax(popt)>numpy.amax(yMean)*0.3:
858 893 # FitGauss=gaus(xSamples,*popt)
859 #
860 # else:
894 #
895 # else:
861 896 # FitGauss=numpy.ones(len(xSamples))*numpy.mean(yMean)
862 897 # print 'Verificador: Dentro', Height
863 898 # except RuntimeError:
864 #
899 #
865 900 # try:
866 901 # for j in range(len(ySamples[1])):
867 902 # yMean2=numpy.append(yMean2,numpy.average([ySamples[1,j],ySamples[2,j]]))
@@ -869,7 +904,7 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
869 904 # FitGauss=gaus(xSamples,*popt)
870 905 # print 'Verificador: Exepcion1', Height
871 906 # except RuntimeError:
872 #
907 #
873 908 # try:
874 909 # popt,pcov = curve_fit(gaus,xSamples,ySamples[1],p0=[1,meanGauss,sigma])
875 910 # FitGauss=gaus(xSamples,*popt)
@@ -880,12 +915,12 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
880 915 # else:
881 916 # FitGauss=numpy.ones(len(xSamples))*numpy.mean(yMean)
882 917 # #print 'Verificador: Fuera', Height
883 #
884 #
885 #
918 #
919 #
920 #
886 921 # Maximun=numpy.amax(yMean)
887 922 # eMinus1=Maximun*numpy.exp(-1)
888 #
923 #
889 924 # HWpos=Find(FitGauss,min(FitGauss, key=lambda value:abs(value-eMinus1)))
890 925 # HalfWidth= xFrec[HWpos]
891 926 # GCpos=Find(FitGauss, numpy.amax(FitGauss))
@@ -894,39 +929,39 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
894 929 # #Vpos=Find(FactNorm, min(FactNorm, key=lambda value:abs(value- numpy.mean(FactNorm) )))
895 930 # #print 'GCpos',GCpos, numpy.amax(FitGauss), 'HWpos',HWpos
896 931 # '''****** Getting Fij ******'''
897 #
932 #
898 933 # GaussCenter=xFrec[GCpos]
899 934 # if (GaussCenter<0 and HalfWidth>0) or (GaussCenter>0 and HalfWidth<0):
900 935 # Fij=abs(GaussCenter)+abs(HalfWidth)+0.0000001
901 936 # else:
902 937 # Fij=abs(GaussCenter-HalfWidth)+0.0000001
903 #
938 #
904 939 # '''****** Getting Frecuency range of significant data ******'''
905 #
940 #
906 941 # Rangpos=Find(FitGauss,min(FitGauss, key=lambda value:abs(value-Maximun*0.10)))
907 #
942 #
908 943 # if Rangpos<GCpos:
909 944 # Range=numpy.array([Rangpos,2*GCpos-Rangpos])
910 945 # else:
911 946 # Range=numpy.array([2*GCpos-Rangpos,Rangpos])
912 #
947 #
913 948 # FrecRange=xFrec[Range[0]:Range[1]]
914 #
949 #
915 950 # #print 'FrecRange', FrecRange
916 951 # '''****** Getting SCPC Slope ******'''
917 #
952 #
918 953 # for i in range(self.nRdPairs):
919 #
954 #
920 955 # if len(FrecRange)>5 and len(FrecRange)<self.nProfiles*0.5:
921 # PhaseRange=moving_average(phase[i,Range[0]:Range[1]],N=3)
922 #
956 # PhaseRange=moving_average(phase[i,Range[0]:Range[1]],N=3)
957 #
923 958 # slope, intercept, r_value, p_value, std_err = stats.linregress(FrecRange,PhaseRange)
924 959 # PhaseSlope[i]=slope
925 960 # PhaseInter[i]=intercept
926 961 # else:
927 962 # PhaseSlope[i]=0
928 963 # PhaseInter[i]=0
929 #
964 #
930 965 # # plt.figure(i+15)
931 966 # # plt.title('FASE ( CH%s*CH%s )' %(self.dataOut.pairsList[i][0],self.dataOut.pairsList[i][1]))
932 967 # # plt.xlabel('Frecuencia (KHz)')
@@ -934,76 +969,76 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
934 969 # # #plt.subplot(311+i)
935 970 # # plt.plot(FrecRange,PhaseRange,'b')
936 971 # # plt.plot(FrecRange,FrecRange*PhaseSlope[i]+PhaseInter[i],'r')
937 #
972 #
938 973 # #plt.axis([-0.6, 0.2, -3.2, 3.2])
939 #
940 #
974 #
975 #
941 976 # '''Getting constant C'''
942 977 # cC=(Fij*numpy.pi)**2
943 #
978 #
944 979 # # '''Getting Eij and Nij'''
945 980 # # (AntennaX0,AntennaY0)=pol2cart(rheader.AntennaCoord0, rheader.AntennaAngl0*numpy.pi/180)
946 981 # # (AntennaX1,AntennaY1)=pol2cart(rheader.AntennaCoord1, rheader.AntennaAngl1*numpy.pi/180)
947 982 # # (AntennaX2,AntennaY2)=pol2cart(rheader.AntennaCoord2, rheader.AntennaAngl2*numpy.pi/180)
948 # #
983 # #
949 984 # # E01=AntennaX0-AntennaX1
950 985 # # N01=AntennaY0-AntennaY1
951 # #
986 # #
952 987 # # E02=AntennaX0-AntennaX2
953 988 # # N02=AntennaY0-AntennaY2
954 # #
989 # #
955 990 # # E12=AntennaX1-AntennaX2
956 991 # # N12=AntennaY1-AntennaY2
957 #
992 #
958 993 # '''****** Getting constants F and G ******'''
959 994 # MijEijNij=numpy.array([[E02,N02], [E12,N12]])
960 995 # MijResult0=(-PhaseSlope[1]*cC) / (2*numpy.pi)
961 # MijResult1=(-PhaseSlope[2]*cC) / (2*numpy.pi)
996 # MijResult1=(-PhaseSlope[2]*cC) / (2*numpy.pi)
962 997 # MijResults=numpy.array([MijResult0,MijResult1])
963 998 # (cF,cG) = numpy.linalg.solve(MijEijNij, MijResults)
964 #
999 #
965 1000 # '''****** Getting constants A, B and H ******'''
966 1001 # W01=numpy.amax(coherence[0])
967 1002 # W02=numpy.amax(coherence[1])
968 1003 # W12=numpy.amax(coherence[2])
969 #
1004 #
970 1005 # WijResult0=((cF*E01+cG*N01)**2)/cC - numpy.log(W01 / numpy.sqrt(numpy.pi/cC))
971 1006 # WijResult1=((cF*E02+cG*N02)**2)/cC - numpy.log(W02 / numpy.sqrt(numpy.pi/cC))
972 1007 # WijResult2=((cF*E12+cG*N12)**2)/cC - numpy.log(W12 / numpy.sqrt(numpy.pi/cC))
973 #
1008 #
974 1009 # WijResults=numpy.array([WijResult0, WijResult1, WijResult2])
975 #
976 # WijEijNij=numpy.array([ [E01**2, N01**2, 2*E01*N01] , [E02**2, N02**2, 2*E02*N02] , [E12**2, N12**2, 2*E12*N12] ])
1010 #
1011 # WijEijNij=numpy.array([ [E01**2, N01**2, 2*E01*N01] , [E02**2, N02**2, 2*E02*N02] , [E12**2, N12**2, 2*E12*N12] ])
977 1012 # (cA,cB,cH) = numpy.linalg.solve(WijEijNij, WijResults)
978 #
1013 #
979 1014 # VxVy=numpy.array([[cA,cH],[cH,cB]])
980 #
1015 #
981 1016 # VxVyResults=numpy.array([-cF,-cG])
982 1017 # (Vx,Vy) = numpy.linalg.solve(VxVy, VxVyResults)
983 1018 # Vzon = Vy
984 1019 # Vmer = Vx
985 1020 # Vmag=numpy.sqrt(Vzon**2+Vmer**2)
986 1021 # Vang=numpy.arctan2(Vmer,Vzon)
987 #
1022 #
988 1023 # if abs(Vy)<100 and abs(Vy)> 0.:
989 1024 # self.dataOut.velocityX=numpy.append(self.dataOut.velocityX, Vzon) #Vmag
990 1025 # #print 'Vmag',Vmag
991 1026 # else:
992 1027 # self.dataOut.velocityX=numpy.append(self.dataOut.velocityX, NaN)
993 #
1028 #
994 1029 # if abs(Vx)<100 and abs(Vx) > 0.:
995 1030 # self.dataOut.velocityY=numpy.append(self.dataOut.velocityY, Vmer) #Vang
996 # #print 'Vang',Vang
1031 # #print 'Vang',Vang
997 1032 # else:
998 1033 # self.dataOut.velocityY=numpy.append(self.dataOut.velocityY, NaN)
999 #
1034 #
1000 1035 # if abs(GaussCenter)<2:
1001 1036 # self.dataOut.velocityV=numpy.append(self.dataOut.velocityV, xFrec[Vpos])
1002 #
1037 #
1003 1038 # else:
1004 1039 # self.dataOut.velocityV=numpy.append(self.dataOut.velocityV, NaN)
1005 #
1006 #
1040 #
1041 #
1007 1042 # # print '********************************************'
1008 1043 # # print 'HalfWidth ', HalfWidth
1009 1044 # # print 'Maximun ', Maximun
@@ -1033,25 +1068,25 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1033 1068 # # print 'PhaseSlope ',PhaseSlope[2]
1034 1069 # # print '********************************************'
1035 1070 # #print 'data_output',shape(self.dataOut.velocityX), shape(self.dataOut.velocityY)
1036 #
1071 #
1037 1072 # #print 'self.dataOut.velocityX', len(self.dataOut.velocityX)
1038 1073 # #print 'self.dataOut.velocityY', len(self.dataOut.velocityY)
1039 1074 # #print 'self.dataOut.velocityV', self.dataOut.velocityV
1040 #
1075 #
1041 1076 # self.data_output[0]=numpy.array(self.dataOut.velocityX)
1042 1077 # self.data_output[1]=numpy.array(self.dataOut.velocityY)
1043 1078 # self.data_output[2]=numpy.array(self.dataOut.velocityV)
1044 #
1079 #
1045 1080 # prin= self.data_output[0][~numpy.isnan(self.data_output[0])]
1046 1081 # print ' '
1047 # print 'VmagAverage',numpy.mean(prin)
1082 # print 'VmagAverage',numpy.mean(prin)
1048 1083 # print ' '
1049 1084 # # plt.figure(5)
1050 1085 # # plt.subplot(211)
1051 1086 # # plt.plot(self.dataOut.velocityX,'yo:')
1052 1087 # # plt.subplot(212)
1053 1088 # # plt.plot(self.dataOut.velocityY,'yo:')
1054 #
1089 #
1055 1090 # # plt.figure(1)
1056 1091 # # # plt.subplot(121)
1057 1092 # # # plt.plot(xFrec,ySamples[0],'k',label='Ch0')
@@ -1060,7 +1095,7 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1060 1095 # # # plt.plot(xFrec,FitGauss,'yo:',label='fit')
1061 1096 # # # plt.legend()
1062 1097 # # plt.title('DATOS A ALTURA DE 2850 METROS')
1063 # #
1098 # #
1064 1099 # # plt.xlabel('Frecuencia (KHz)')
1065 1100 # # plt.ylabel('Magnitud')
1066 1101 # # # plt.subplot(122)
@@ -1070,7 +1105,7 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1070 1105 # # plt.plot(xFrec,FactNorm)
1071 1106 # # plt.axis([-4, 4, 0, 0.15])
1072 1107 # # # plt.xlabel('SelfSpectra KHz')
1073 # #
1108 # #
1074 1109 # # plt.figure(10)
1075 1110 # # # plt.subplot(121)
1076 1111 # # plt.plot(xFrec,ySamples[0],'b',label='Ch0')
@@ -1079,7 +1114,7 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1079 1114 # # # plt.plot(xFrec,FitGauss,'yo:',label='fit')
1080 1115 # # plt.legend()
1081 1116 # # plt.title('SELFSPECTRA EN CANALES')
1082 # #
1117 # #
1083 1118 # # plt.xlabel('Frecuencia (KHz)')
1084 1119 # # plt.ylabel('Magnitud')
1085 1120 # # # plt.subplot(122)
@@ -1089,19 +1124,19 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1089 1124 # # # plt.plot(xFrec,FactNorm)
1090 1125 # # # plt.axis([-4, 4, 0, 0.15])
1091 1126 # # # plt.xlabel('SelfSpectra KHz')
1092 # #
1127 # #
1093 1128 # # plt.figure(9)
1094 # #
1095 # #
1129 # #
1130 # #
1096 1131 # # plt.title('DATOS SUAVIZADOS')
1097 1132 # # plt.xlabel('Frecuencia (KHz)')
1098 1133 # # plt.ylabel('Magnitud')
1099 1134 # # plt.plot(xFrec,SmoothSPC,'g')
1100 # #
1135 # #
1101 1136 # # #plt.plot(xFrec,FactNorm)
1102 1137 # # plt.axis([-4, 4, 0, 0.15])
1103 1138 # # # plt.xlabel('SelfSpectra KHz')
1104 # # #
1139 # # #
1105 1140 # # plt.figure(2)
1106 1141 # # # #plt.subplot(121)
1107 1142 # # plt.plot(xFrec,yMean,'r',label='Mean SelfSpectra')
@@ -1115,7 +1150,7 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1115 1150 # # plt.xlabel('Frecuencia (KHz)')
1116 1151 # # plt.ylabel('Magnitud')
1117 1152 # # plt.legend()
1118 # # #
1153 # # #
1119 1154 # # # plt.figure(3)
1120 1155 # # # plt.subplot(311)
1121 1156 # # # #plt.plot(xFrec,phase[0])
@@ -1125,30 +1160,23 class BLTRSpectraReader (ProcessingUnit, FileHeaderBLTR, RecordHeaderBLTR, JRODa
1125 1160 # # # plt.subplot(313)
1126 1161 # # # plt.plot(xFrec,phase[2],'g')
1127 1162 # # # #plt.plot(xFrec,phase[2])
1128 # # #
1163 # # #
1129 1164 # # # plt.figure(4)
1130 # # #
1165 # # #
1131 1166 # # # plt.plot(xSamples,coherence[0],'b')
1132 1167 # # # plt.plot(xSamples,coherence[1],'r')
1133 1168 # # # plt.plot(xSamples,coherence[2],'g')
1134 1169 # # plt.show()
1135 # # #
1170 # # #
1136 1171 # # # plt.clf()
1137 1172 # # # plt.cla()
1138 # # # plt.close()
1139 #
1140 # print ' '
1141
1142
1143
1144 self.BlockCounter+=2
1145
1173 # # # plt.close()
1174 #
1175 # print ' '
1176
1177 self.BlockCounter += 2
1178
1146 1179 else:
1147 self.fileSelector+=1
1148 self.BlockCounter=0
1180 self.fileSelector += 1
1181 self.BlockCounter = 0
1149 1182 print "Next File"
1150
1151
1152
1153
1154
Show file before | Show file after | Hide comments
This diff has been collapsed as it changes many lines, (1027 lines changed) Show them Hide them
@@ -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,11 +7,9 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 from scipy import asarray as ar,exp
12 from scipy import asarray as ar, exp
14 13 from scipy import stats
15 14
16 15 from numpy.ma.core import getdata
@@ -30,152 +29,168 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
39
36
40 37 def read(self):
41
38
42 39 raise NotImplementedError
43
40
44 41 def write(self):
45
42
46 43 raise NotImplementedError
47
44
48 45 def printInfo(self):
49
50 message = "#"*50 + "\n"
46
47 message = "#" * 50 + "\n"
51 48 message += self.__class__.__name__.upper() + "\n"
52 message += "#"*50 + "\n"
53
49 message += "#" * 50 + "\n"
50
54 51 keyList = self.__dict__.keys()
55 52 keyList.sort()
56
53
57 54 for key in keyList:
58 message += "%s = %s" %(key, self.__dict__[key]) + "\n"
59
55 message += "%s = %s" % (key, self.__dict__[key]) + "\n"
56
60 57 if "size" not in keyList:
61 58 attr = getattr(self, "size")
62
59
63 60 if attr:
64 message += "%s = %s" %("size", attr) + "\n"
65
66 #print message
67
68
69 FILE_HEADER = numpy.dtype([ #HEADER 1024bytes
70 ('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
75 ('Hdummy',numpy.str_,512), #Reserved space
76 #Main chunk 8bytes
77 ('Msign',numpy.str_,4), #Main chunk signature FZKF or NUIG
78 ('MsizeData','<i4'), #Size of data block main chunk
79 #Processing DSP parameters 36bytes
80 ('PPARsign',numpy.str_,4), #PPAR signature
81 ('PPARsize','<i4'), #PPAR size of block
82 ('PPARprf','<i4'), #Pulse repetition frequency
83 ('PPARpdr','<i4'), #Pulse duration
84 ('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.
89 #Service DSP parameters 112bytes
90 ('SPARatt','<i4'), #STC attenuation on the lowest ranges on/off
91 ('SPARtx','<i4'), #OBSOLETE
92 ('SPARaddGain0','<f4'), #OBSOLETE
93 ('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
98 ('SPARcal','<i4'), #OBSOLETE
99 ('SPARnos','<i4'), #OBSOLETE
100 ('SPARof0','<i4'), #detection threshold
101 ('SPARof1','<i4'), #OBSOLETE
102 ('SPARswt','<i4'), #2nd moment estimation threshold
103 ('SPARsum','<i4'), #OBSOLETE
104 ('SPARosc','<i4'), #flag Oscillosgram mode
105 ('SPARtst','<i4'), #OBSOLETE
106 ('SPARcor','<i4'), #OBSOLETE
107 ('SPARofs','<i4'), #OBSOLETE
108 ('SPARhsn','<i4'), #Hildebrand div noise detection on noise gate
109 ('SPARhsa','<f4'), #Hildebrand div noise detection on all gates
110 ('SPARcalibPow_M','<f4'), #OBSOLETE
111 ('SPARcalibSNR_M','<f4'), #OBSOLETE
112 ('SPARcalibPow_S','<f4'), #OBSOLETE
113 ('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
117 ('SPARprc','<i4'),]) #flag - Moment estimation switched on/off
118
119
120
61 message += "%s = %s" % ("size", attr) + "\n"
62
63 # print message
64
65
66 FILE_HEADER = numpy.dtype([ # HEADER 1024bytes
67 ('Hname', 'a32'), # Original file name
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),
76 ('Hdummy', numpy.str_, 512), # Reserved space
77 # Main chunk 8bytes
78 # Main chunk signature FZKF or NUIG
79 ('Msign', numpy.str_, 4),
80 ('MsizeData', '<i4'), # Size of data block main chunk
81 # Processing DSP parameters 36bytes
82 ('PPARsign', numpy.str_, 4), # PPAR signature
83 ('PPARsize', '<i4'), # PPAR size of block
84 ('PPARprf', '<i4'), # Pulse repetition frequency
85 ('PPARpdr', '<i4'), # Pulse duration
86 ('PPARsft', '<i4'), # FFT length
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'),
95 # Service DSP parameters 112bytes
96 # STC attenuation on the lowest ranges on/off
97 ('SPARatt', '<i4'),
98 ('SPARtx', '<i4'), # OBSOLETE
99 ('SPARaddGain0', '<f4'), # OBSOLETE
100 ('SPARaddGain1', '<f4'), # 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'),
109 ('SPARcal', '<i4'), # OBSOLETE
110 ('SPARnos', '<i4'), # OBSOLETE
111 ('SPARof0', '<i4'), # detection threshold
112 ('SPARof1', '<i4'), # OBSOLETE
113 ('SPARswt', '<i4'), # 2nd moment estimation threshold
114 ('SPARsum', '<i4'), # OBSOLETE
115 ('SPARosc', '<i4'), # flag Oscillosgram mode
116 ('SPARtst', '<i4'), # OBSOLETE
117 ('SPARcor', '<i4'), # OBSOLETE
118 ('SPARofs', '<i4'), # OBSOLETE
119 # Hildebrand div noise detection on noise gate
120 ('SPARhsn', '<i4'),
121 # Hildebrand div noise detection on all gates
122 ('SPARhsa', '<f4'),
123 ('SPARcalibPow_M', '<f4'), # OBSOLETE
124 ('SPARcalibSNR_M', '<f4'), # OBSOLETE
125 ('SPARcalibPow_S', '<f4'), # OBSOLETE
126 ('SPARcalibSNR_S', '<f4'), # OBSOLETE
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'),
133 ('SPARprc', '<i4'), ]) # flag - Moment estimation switched on/off
134
135
121 136 class FileHeaderMIRA35c(Header):
122
137
123 138 def __init__(self):
124
125 self.Hname= None
126 self.Htime= None
127 self.Hoper= None
128 self.Hplace= None
129 self.Hdescr= None
130 self.Hdummy= None
131
132 self.Msign=None
133 self.MsizeData=None
134
135 self.PPARsign=None
136 self.PPARsize=None
137 self.PPARprf=None
138 self.PPARpdr=None
139 self.PPARsft=None
140 self.PPARavc=None
141 self.PPARihp=None
142 self.PPARchg=None
143 self.PPARpol=None
144 #Service DSP parameters
145 self.SPARatt=None
146 self.SPARtx=None
147 self.SPARaddGain0=None
148 self.SPARaddGain1=None
149 self.SPARwnd=None
150 self.SPARpos=None
151 self.SPARadd=None
152 self.SPARlen=None
153 self.SPARcal=None
154 self.SPARnos=None
155 self.SPARof0=None
156 self.SPARof1=None
157 self.SPARswt=None
158 self.SPARsum=None
159 self.SPARosc=None
160 self.SPARtst=None
161 self.SPARcor=None
162 self.SPARofs=None
163 self.SPARhsn=None
164 self.SPARhsa=None
165 self.SPARcalibPow_M=None
166 self.SPARcalibSNR_M=None
167 self.SPARcalibPow_S=None
168 self.SPARcalibSNR_S=None
169 self.SPARrawGate1=None
170 self.SPARrawGate2=None
171 self.SPARraw=None
172 self.SPARprc=None
173
174 self.FHsize=1180
175
139
140 self.Hname = None
141 self.Htime = None
142 self.Hoper = None
143 self.Hplace = None
144 self.Hdescr = None
145 self.Hdummy = None
146
147 self.Msign = None
148 self.MsizeData = None
149
150 self.PPARsign = None
151 self.PPARsize = None
152 self.PPARprf = None
153 self.PPARpdr = None
154 self.PPARsft = None
155 self.PPARavc = None
156 self.PPARihp = None
157 self.PPARchg = None
158 self.PPARpol = None
159 # Service DSP parameters
160 self.SPARatt = None
161 self.SPARtx = None
162 self.SPARaddGain0 = None
163 self.SPARaddGain1 = None
164 self.SPARwnd = None
165 self.SPARpos = None
166 self.SPARadd = None
167 self.SPARlen = None
168 self.SPARcal = None
169 self.SPARnos = None
170 self.SPARof0 = None
171 self.SPARof1 = None
172 self.SPARswt = None
173 self.SPARsum = None
174 self.SPARosc = None
175 self.SPARtst = None
176 self.SPARcor = None
177 self.SPARofs = None
178 self.SPARhsn = None
179 self.SPARhsa = None
180 self.SPARcalibPow_M = None
181 self.SPARcalibSNR_M = None
182 self.SPARcalibPow_S = None
183 self.SPARcalibSNR_S = None
184 self.SPARrawGate1 = None
185 self.SPARrawGate2 = None
186 self.SPARraw = None
187 self.SPARprc = None
188
189 self.FHsize = 1180
190
176 191 def FHread(self, fp):
177
178 header = numpy.fromfile(fp, FILE_HEADER,1)
192
193 header = numpy.fromfile(fp, FILE_HEADER, 1)
179 194 ''' numpy.fromfile(file, dtype, count, sep='')
180 195 file : file or str
181 196 Open file object or filename.
@@ -194,85 +209,83 class FileHeaderMIRA35c(Header):
194 209 at least one whitespace.
195 210
196 211 '''
197
198
199 self.Hname= str(header['Hname'][0])
200 self.Htime= str(header['Htime'][0])
201 self.Hoper= str(header['Hoper'][0])
202 self.Hplace= str(header['Hplace'][0])
203 self.Hdescr= str(header['Hdescr'][0])
204 self.Hdummy= str(header['Hdummy'][0])
205 #1024
206
207 self.Msign=str(header['Msign'][0])
208 self.MsizeData=header['MsizeData'][0]
209 #8
210
211 self.PPARsign=str(header['PPARsign'][0])
212 self.PPARsize=header['PPARsize'][0]
213 self.PPARprf=header['PPARprf'][0]
214 self.PPARpdr=header['PPARpdr'][0]
215 self.PPARsft=header['PPARsft'][0]
216 self.PPARavc=header['PPARavc'][0]
217 self.PPARihp=header['PPARihp'][0]
218 self.PPARchg=header['PPARchg'][0]
219 self.PPARpol=header['PPARpol'][0]
220 #Service DSP parameters
221 #36
222
223 self.SPARatt=header['SPARatt'][0]
224 self.SPARtx=header['SPARtx'][0]
225 self.SPARaddGain0=header['SPARaddGain0'][0]
226 self.SPARaddGain1=header['SPARaddGain1'][0]
227 self.SPARwnd=header['SPARwnd'][0]
228 self.SPARpos=header['SPARpos'][0]
229 self.SPARadd=header['SPARadd'][0]
230 self.SPARlen=header['SPARlen'][0]
231 self.SPARcal=header['SPARcal'][0]
232 self.SPARnos=header['SPARnos'][0]
233 self.SPARof0=header['SPARof0'][0]
234 self.SPARof1=header['SPARof1'][0]
235 self.SPARswt=header['SPARswt'][0]
236 self.SPARsum=header['SPARsum'][0]
237 self.SPARosc=header['SPARosc'][0]
238 self.SPARtst=header['SPARtst'][0]
239 self.SPARcor=header['SPARcor'][0]
240 self.SPARofs=header['SPARofs'][0]
241 self.SPARhsn=header['SPARhsn'][0]
242 self.SPARhsa=header['SPARhsa'][0]
243 self.SPARcalibPow_M=header['SPARcalibPow_M'][0]
244 self.SPARcalibSNR_M=header['SPARcalibSNR_M'][0]
245 self.SPARcalibPow_S=header['SPARcalibPow_S'][0]
246 self.SPARcalibSNR_S=header['SPARcalibSNR_S'][0]
247 self.SPARrawGate1=header['SPARrawGate1'][0]
248 self.SPARrawGate2=header['SPARrawGate2'][0]
249 self.SPARraw=header['SPARraw'][0]
250 self.SPARprc=header['SPARprc'][0]
251 #112
252 #1180
253 #print 'Pointer fp header', fp.tell()
254 #print ' '
255 #print 'SPARrawGate'
256 #print self.SPARrawGate2 - self.SPARrawGate1
257
258 #print ' '
259 #print 'Hname'
260 #print self.Hname
261
262 #print ' '
263 #print 'Msign'
264 #print self.Msign
265
212
213 self.Hname = str(header['Hname'][0])
214 self.Htime = str(header['Htime'][0])
215 self.Hoper = str(header['Hoper'][0])
216 self.Hplace = str(header['Hplace'][0])
217 self.Hdescr = str(header['Hdescr'][0])
218 self.Hdummy = str(header['Hdummy'][0])
219 # 1024
220
221 self.Msign = str(header['Msign'][0])
222 self.MsizeData = header['MsizeData'][0]
223 # 8
224
225 self.PPARsign = str(header['PPARsign'][0])
226 self.PPARsize = header['PPARsize'][0]
227 self.PPARprf = header['PPARprf'][0]
228 self.PPARpdr = header['PPARpdr'][0]
229 self.PPARsft = header['PPARsft'][0]
230 self.PPARavc = header['PPARavc'][0]
231 self.PPARihp = header['PPARihp'][0]
232 self.PPARchg = header['PPARchg'][0]
233 self.PPARpol = header['PPARpol'][0]
234 # Service DSP parameters
235 # 36
236
237 self.SPARatt = header['SPARatt'][0]
238 self.SPARtx = header['SPARtx'][0]
239 self.SPARaddGain0 = header['SPARaddGain0'][0]
240 self.SPARaddGain1 = header['SPARaddGain1'][0]
241 self.SPARwnd = header['SPARwnd'][0]
242 self.SPARpos = header['SPARpos'][0]
243 self.SPARadd = header['SPARadd'][0]
244 self.SPARlen = header['SPARlen'][0]
245 self.SPARcal = header['SPARcal'][0]
246 self.SPARnos = header['SPARnos'][0]
247 self.SPARof0 = header['SPARof0'][0]
248 self.SPARof1 = header['SPARof1'][0]
249 self.SPARswt = header['SPARswt'][0]
250 self.SPARsum = header['SPARsum'][0]
251 self.SPARosc = header['SPARosc'][0]
252 self.SPARtst = header['SPARtst'][0]
253 self.SPARcor = header['SPARcor'][0]
254 self.SPARofs = header['SPARofs'][0]
255 self.SPARhsn = header['SPARhsn'][0]
256 self.SPARhsa = header['SPARhsa'][0]
257 self.SPARcalibPow_M = header['SPARcalibPow_M'][0]
258 self.SPARcalibSNR_M = header['SPARcalibSNR_M'][0]
259 self.SPARcalibPow_S = header['SPARcalibPow_S'][0]
260 self.SPARcalibSNR_S = header['SPARcalibSNR_S'][0]
261 self.SPARrawGate1 = header['SPARrawGate1'][0]
262 self.SPARrawGate2 = header['SPARrawGate2'][0]
263 self.SPARraw = header['SPARraw'][0]
264 self.SPARprc = header['SPARprc'][0]
265 # 112
266 # 1180
267 # print 'Pointer fp header', fp.tell()
268 # print ' '
269 # print 'SPARrawGate'
270 # print self.SPARrawGate2 - self.SPARrawGate1
271
272 # print ' '
273 # print 'Hname'
274 # print self.Hname
275
276 # print ' '
277 # print 'Msign'
278 # print self.Msign
279
266 280 def write(self, fp):
267
281
268 282 headerTuple = (self.Hname,
269 283 self.Htime,
270 284 self.Hoper,
271 285 self.Hplace,
272 286 self.Hdescr,
273 287 self.Hdummy)
274
275
288
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)
@@ -290,72 +303,70 class FileHeaderMIRA35c(Header):
290 303 first converting it to the closest Python type, and then using "format" % item.
291 304
292 305 '''
293
306
294 307 return 1
295 308
296 SRVI_HEADER = numpy.dtype([
297 ('SignatureSRVI1',numpy.str_,4),#
298 ('SizeOfDataBlock1','<i4'),#
299 ('DataBlockTitleSRVI1',numpy.str_,4),#
300 ('SizeOfSRVI1','<i4'),])#
301
302 class SRVIHeader(Header):
309
310 SRVI_HEADER = numpy.dtype([
311 ('SignatureSRVI1', numpy.str_, 4),
312 ('SizeOfDataBlock1', '<i4'),
313 ('DataBlockTitleSRVI1', numpy.str_, 4),
314 ('SizeOfSRVI1', '<i4'), ])
315
316
317 class SRVIHeader(Header):
303 318 def __init__(self, SignatureSRVI1=0, SizeOfDataBlock1=0, DataBlockTitleSRVI1=0, SizeOfSRVI1=0):
304
319
305 320 self.SignatureSRVI1 = SignatureSRVI1
306 321 self.SizeOfDataBlock1 = SizeOfDataBlock1
307 322 self.DataBlockTitleSRVI1 = DataBlockTitleSRVI1
308 self.SizeOfSRVI1 = SizeOfSRVI1
309
310 self.SRVIHsize=16
311
312 def SRVIread(self, fp):
313
314 header = numpy.fromfile(fp, SRVI_HEADER,1)
323 self.SizeOfSRVI1 = SizeOfSRVI1
324
325 self.SRVIHsize = 16
326
327 def SRVIread(self, fp):
328
329 header = numpy.fromfile(fp, SRVI_HEADER, 1)
315 330
316 331 self.SignatureSRVI1 = str(header['SignatureSRVI1'][0])
317 332 self.SizeOfDataBlock1 = header['SizeOfDataBlock1'][0]
318 333 self.DataBlockTitleSRVI1 = str(header['DataBlockTitleSRVI1'][0])
319 334 self.SizeOfSRVI1 = header['SizeOfSRVI1'][0]
320 #16
335 # 16
321 336 print 'Pointer fp SRVIheader', fp.tell()
322
323
324 SRVI_STRUCTURE = numpy.dtype([
325 ('frame_cnt','<u4'),#
326 ('time_t','<u4'), #
327 ('tpow','<f4'), #
328 ('npw1','<f4'), #
329 ('npw2','<f4'), #
330 ('cpw1','<f4'), #
331 ('pcw2','<f4'), #
332 ('ps_err','<u4'), #
333 ('te_err','<u4'), #
334 ('rc_err','<u4'), #
335 ('grs1','<u4'), #
336 ('grs2','<u4'), #
337 ('azipos','<f4'), #
338 ('azivel','<f4'), #
339 ('elvpos','<f4'), #
340 ('elvvel','<f4'), #
341 ('northAngle','<f4'), #
342 ('microsec','<u4'), #
343 ('azisetvel','<f4'), #
344 ('elvsetpos','<f4'), #
345 ('RadarConst','<f4'),]) #
346 337
347 338
339 SRVI_STRUCTURE = numpy.dtype([
340 ('frame_cnt', '<u4'),
341 ('time_t', '<u4'), #
342 ('tpow', '<f4'), #
343 ('npw1', '<f4'), #
344 ('npw2', '<f4'), #
345 ('cpw1', '<f4'), #
346 ('pcw2', '<f4'), #
347 ('ps_err', '<u4'), #
348 ('te_err', '<u4'), #
349 ('rc_err', '<u4'), #
350 ('grs1', '<u4'), #
351 ('grs2', '<u4'), #
352 ('azipos', '<f4'), #
353 ('azivel', '<f4'), #
354 ('elvpos', '<f4'), #
355 ('elvvel', '<f4'), #
356 ('northAngle', '<f4'),
357 ('microsec', '<u4'), #
358 ('azisetvel', '<f4'), #
359 ('elvsetpos', '<f4'), #
360 ('RadarConst', '<f4'), ]) #
348 361
349 362
350 363 class RecordHeader(Header):
351
352
353 def __init__(self, frame_cnt=0, time_t= 0, tpow=0, npw1=0, npw2=0,
354 cpw1=0, pcw2=0, ps_err=0, te_err=0, rc_err=0, grs1=0,
355 grs2=0, azipos=0, azivel=0, elvpos=0, elvvel=0, northangle=0,
356 microsec=0, azisetvel=0, elvsetpos=0, RadarConst=0 , RecCounter=0, Off2StartNxtRec=0):
357
358
364
365 def __init__(self, frame_cnt=0, time_t=0, tpow=0, npw1=0, npw2=0,
366 cpw1=0, pcw2=0, ps_err=0, te_err=0, rc_err=0, grs1=0,
367 grs2=0, azipos=0, azivel=0, elvpos=0, elvvel=0, northangle=0,
368 microsec=0, azisetvel=0, elvsetpos=0, RadarConst=0, RecCounter=0, Off2StartNxtRec=0):
369
359 370 self.frame_cnt = frame_cnt
360 371 self.dwell = time_t
361 372 self.tpow = tpow
@@ -377,22 +388,22 class RecordHeader(Header):
377 388 self.azisetvel = azisetvel
378 389 self.elvsetpos = elvsetpos
379 390 self.RadarConst = RadarConst
380 self.RHsize=84
391 self.RHsize = 84
381 392 self.RecCounter = RecCounter
382 self.Off2StartNxtRec=Off2StartNxtRec
383
393 self.Off2StartNxtRec = Off2StartNxtRec
394
384 395 def RHread(self, fp):
385
386 #startFp = open(fp,"rb") #The method tell() returns the current position of the file read/write pointer within the file.
387
396
397 # startFp = open(fp,"rb") #The method tell() returns the current position of the file read/write pointer within the file.
398
388 399 #OffRHeader= 1180 + self.RecCounter*(self.Off2StartNxtRec)
389 400 #startFp.seek(OffRHeader, os.SEEK_SET)
390
391 #print 'Posicion del bloque: ',OffRHeader
392
393 header = numpy.fromfile(fp,SRVI_STRUCTURE,1)
394
395 self.frame_cnt = header['frame_cnt'][0]#
401
402 # print 'Posicion del bloque: ',OffRHeader
403
404 header = numpy.fromfile(fp, SRVI_STRUCTURE, 1)
405
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] #
@@ -413,26 +424,26 class RecordHeader(Header):
413 424 self.azisetvel = header['azisetvel'][0] #
414 425 self.elvsetpos = header['elvsetpos'][0] #
415 426 self.RadarConst = header['RadarConst'][0] #
416 #84
417
418 #print 'Pointer fp RECheader', fp.tell()
419
427 # 84
428
429 # print 'Pointer fp RECheader', fp.tell()
430
420 431 #self.ipp= 0.5*(SPEED_OF_LIGHT/self.PRFhz)
421
432
422 433 #self.RHsize = 180+20*self.nChannels
423 434 #self.Datasize= self.nProfiles*self.nChannels*self.nHeights*2*4
424 #print 'Datasize',self.Datasize
435 # print 'Datasize',self.Datasize
425 436 #endFp = self.OffsetStartHeader + self.RecCounter*self.Off2StartNxtRec
426
437
427 438 print '=============================================='
428
439
429 440 print '=============================================='
430 441
431
432 442 return 1
433
434 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
435
443
444
445 class MIRA35CReader (ProcessingUnit, FileHeaderMIRA35c, SRVIHeader, RecordHeader):
446
436 447 path = None
437 448 startDate = None
438 449 endDate = None
@@ -440,87 +451,84 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
440 451 endTime = None
441 452 walk = None
442 453 isConfig = False
443
444
445 fileList= None
446
447 #metadata
448 TimeZone= None
449 Interval= None
450 heightList= None
451
452 #data
453 data= None
454 utctime= None
455
456
457
458 def __init__(self, **kwargs):
459
460 #Eliminar de la base la herencia
454
455 fileList = None
456
457 # metadata
458 TimeZone = None
459 Interval = None
460 heightList = None
461
462 # data
463 data = None
464 utctime = None
465
466 def __init__(self, **kwargs):
467
468 # Eliminar de la base la herencia
461 469 ProcessingUnit.__init__(self, **kwargs)
462 470 self.PointerReader = 0
463 471 self.FileHeaderFlag = False
464 472 self.utc = None
465 473 self.ext = ".zspca"
466 474 self.optchar = "P"
467 self.fpFile=None
475 self.fpFile = None
468 476 self.fp = None
469 self.BlockCounter=0
477 self.BlockCounter = 0
470 478 self.dtype = None
471 479 self.fileSizeByHeader = None
472 480 self.filenameList = []
473 481 self.fileSelector = 0
474 self.Off2StartNxtRec=0
475 self.RecCounter=0
482 self.Off2StartNxtRec = 0
483 self.RecCounter = 0
476 484 self.flagNoMoreFiles = 0
477 self.data_spc=None
478 #self.data_cspc=None
479 self.data_output=None
485 self.data_spc = None
486 # self.data_cspc=None
487 self.data_output = None
480 488 self.path = None
481 self.OffsetStartHeader=0
482 self.Off2StartData=0
489 self.OffsetStartHeader = 0
490 self.Off2StartData = 0
483 491 self.ipp = 0
484 self.nFDTdataRecors=0
492 self.nFDTdataRecors = 0
485 493 self.blocksize = 0
486 494 self.dataOut = Spectra()
487 self.profileIndex = 1 #Always
488 self.dataOut.flagNoData=False
495 self.profileIndex = 1 # Always
496 self.dataOut.flagNoData = False
489 497 self.dataOut.nRdPairs = 0
490 498 self.dataOut.pairsList = []
491 self.dataOut.data_spc=None
492
493 self.dataOut.normFactor=1
499 self.dataOut.data_spc = None
500
501 self.dataOut.normFactor = 1
494 502 self.nextfileflag = True
495 503 self.dataOut.RadarConst = 0
496 504 self.dataOut.HSDV = []
497 505 self.dataOut.NPW = []
498 506 self.dataOut.COFA = []
499 507 self.dataOut.noise = 0
500
501
508
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.
505 512 It also creates an organized list with the names of the files to read.
506 513 '''
507 #self.__checkPath()
508
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
511 nFiles=0 #File Counter
512 FileList=[] #A list is created that will contain the .fdt files
513 for IndexFile in ListaData :
514 if '.zspca' in IndexFile and '.gz' not in IndexFile:
514 # self.__checkPath()
515
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)
520 nFiles = 0 # File Counter
521 FileList = [] # A list is created that will contain the .fdt files
522 for IndexFile in ListaData:
523 if '.zspca' in IndexFile and '.gz' not in IndexFile:
515 524 FileList.append(IndexFile)
516 nFiles+=1
517
518 #print 'Files2Read'
519 #print 'Existen '+str(nFiles)+' archivos .fdt'
520
521 self.filenameList=FileList #List of files from least to largest by names
522
523
525 nFiles += 1
526
527 # print 'Files2Read'
528 # print 'Existen '+str(nFiles)+' archivos .fdt'
529
530 self.filenameList = FileList # List of files from least to largest by names
531
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.
@@ -530,274 +538,265 class MIRA35CReader (ProcessingUnit,FileHeaderMIRA35c,SRVIHeader,RecordHeader):
530 538 if not self.isConfig:
531 539 self.setup(**kwargs)
532 540 self.isConfig = True
533
541
534 542 self.getData()
535
536
543
537 544 def setup(self, path=None,
538 startDate=None,
539 endDate=None,
540 startTime=None,
541 endTime=None,
542 walk=True,
543 timezone='utc',
544 code = None,
545 online=False,
546 ReadMode=None, **kwargs):
547
545 startDate=None,
546 endDate=None,
547 startTime=None,
548 endTime=None,
549 walk=True,
550 timezone='utc',
551 code=None,
552 online=False,
553 ReadMode=None, **kwargs):
554
548 555 self.isConfig = True
549
550 self.path=path
551 self.startDate=startDate
552 self.endDate=endDate
553 self.startTime=startTime
554 self.endTime=endTime
555 self.walk=walk
556 #self.ReadMode=int(ReadMode)
557
556
557 self.path = path
558 self.startDate = startDate
559 self.endDate = endDate
560 self.startTime = startTime
561 self.endTime = endTime
562 self.walk = walk
563 # self.ReadMode=int(ReadMode)
564
558 565 pass
559
560
566
561 567 def getData(self):
562 568 '''
563 569 Before starting this function, you should check that there is still an unread file,
564 570 If there are still blocks to read or if the data block is empty.
565
571
566 572 You should call the file "read".
567
573
568 574 '''
569
575
570 576 if self.flagNoMoreFiles:
571 577 self.dataOut.flagNoData = True
572 578 print 'NoData se vuelve true'
573 579 return 0
574
575 self.fp=self.path
580
581 self.fp = self.path
576 582 self.Files2Read(self.fp)
577 583 self.readFile(self.fp)
578
579 self.dataOut.data_spc = self.dataOut_spc#self.data_spc.copy()
584
585 self.dataOut.data_spc = self.dataOut_spc # self.data_spc.copy()
580 586 self.dataOut.RadarConst = self.RadarConst
581 self.dataOut.data_output=self.data_output
587 self.dataOut.data_output = self.data_output
582 588 self.dataOut.noise = self.dataOut.getNoise()
583 #print 'ACAAAAAA', self.dataOut.noise
584 self.dataOut.data_spc = self.dataOut.data_spc+self.dataOut.noise
585 #print 'self.dataOut.noise',self.dataOut.noise
586
587
589 # print 'ACAAAAAA', self.dataOut.noise
590 self.dataOut.data_spc = self.dataOut.data_spc + self.dataOut.noise
591 # print 'self.dataOut.noise',self.dataOut.noise
592
588 593 return self.dataOut.data_spc
589
590
591 def readFile(self,fp):
594
595 def readFile(self, fp):
592 596 '''
593 597 You must indicate if you are reading in Online or Offline mode and load the
594 598 The parameters for this file reading mode.
595
599
596 600 Then you must do 2 actions:
597
601
598 602 1. Get the BLTR FileHeader.
599 603 2. Start reading the first block.
600 604 '''
601
602 #The address of the folder is generated the name of the .fdt file that will be read
603 print "File: ",self.fileSelector+1
604
605
606 # The address of the folder is generated the name of the .fdt file that will be read
607 print "File: ", self.fileSelector + 1
608
605 609 if self.fileSelector < len(self.filenameList):
606
607 self.fpFile=str(fp)+'/'+str(self.filenameList[self.fileSelector])
608
609 if self.nextfileflag==True:
610 self.fp = open(self.fpFile,"rb")
611 self.nextfileflag==False
612
610
611 self.fpFile = str(fp) + '/' + \
612 str(self.filenameList[self.fileSelector])
613
614 if self.nextfileflag == True:
615 self.fp = open(self.fpFile, "rb")
616 self.nextfileflag == False
617
613 618 '''HERE STARTING THE FILE READING'''
614
615
619
616 620 self.fheader = FileHeaderMIRA35c()
617 self.fheader.FHread(self.fp) #Bltr FileHeader Reading
618
621 self.fheader.FHread(self.fp) # Bltr FileHeader Reading
619 622
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
626
628
627 629 self.Num_inCoh = self.fheader.PPARavc
628 630 self.dataOut.PRF = self.fheader.PPARprf
629 self.dataOut.frequency = 34.85*10**9
630 self.Lambda = SPEED_OF_LIGHT/self.dataOut.frequency
631 self.dataOut.ippSeconds= 1./float(self.dataOut.PRF)
632
633 pulse_width = self.fheader.PPARpdr * 10**-9
631 self.dataOut.frequency = 34.85 * 10**9
632 self.Lambda = SPEED_OF_LIGHT / self.dataOut.frequency
633 self.dataOut.ippSeconds = 1. / float(self.dataOut.PRF)
634
635 pulse_width = self.fheader.PPARpdr * 10**-9
634 636 self.__deltaHeigth = 0.5 * SPEED_OF_LIGHT * pulse_width
635
636 self.data_spc = numpy.zeros((self.Num_Hei, self.Num_Bins,2))#
637
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 self.ETA = numpy.zeros(([2,self.Num_Hei]))
641
642
643
644 self.readBlock() #Block reading
645
642 self.ETA = numpy.zeros(([2, self.Num_Hei]))
643
644 self.readBlock() # Block reading
645
646 646 else:
647 647 print 'readFile FlagNoData becomes true'
648 self.flagNoMoreFiles=True
648 self.flagNoMoreFiles = True
649 649 self.dataOut.flagNoData = True
650 650 self.FileHeaderFlag == True
651 return 0
652
653
654
655 def readBlock(self):
651 return 0
652
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.
658
656
659 657 Then the following is done:
660
658
661 659 1. Read the RecordHeader
662 660 2. Fill the buffer with the current block number.
663
661
664 662 '''
665
663
666 664 if self.PointerReader > 1180:
667 self.fp.seek(self.PointerReader , os.SEEK_SET)
665 self.fp.seek(self.PointerReader, os.SEEK_SET)
668 666 self.FirstPoint = self.PointerReader
669
670 else :
667
668 else:
671 669 self.FirstPoint = 1180
672
673
674
670
675 671 self.srviHeader = SRVIHeader()
676
677 self.srviHeader.SRVIread(self.fp) #Se obtiene la cabecera del SRVI
678
679 self.blocksize = self.srviHeader.SizeOfDataBlock1 # Se obtiene el tamao del bloque
680
672
673 self.srviHeader.SRVIread(self.fp) # Se obtiene la cabecera del SRVI
674
675 self.blocksize = self.srviHeader.SizeOfDataBlock1 # Se obtiene el tamao del bloque
676
681 677 if self.blocksize == 148:
682 678 print 'blocksize == 148 bug'
683 jump = numpy.fromfile(self.fp,[('jump',numpy.str_,140)] ,1)
684
685 self.srviHeader.SRVIread(self.fp) #Se obtiene la cabecera del SRVI
686
679 jump = numpy.fromfile(self.fp, [('jump', numpy.str_, 140)], 1)
680
681 # Se obtiene la cabecera del SRVI
682 self.srviHeader.SRVIread(self.fp)
683
687 684 if not self.srviHeader.SizeOfSRVI1:
688 self.fileSelector+=1
689 self.nextfileflag==True
685 self.fileSelector += 1
686 self.nextfileflag == True
690 687 self.FileHeaderFlag == True
691
688
692 689 self.recordheader = RecordHeader()
693 690 self.recordheader.RHread(self.fp)
694 691 self.RadarConst = self.recordheader.RadarConst
695 692 dwell = self.recordheader.time_t
696 693 npw1 = self.recordheader.npw1
697 694 npw2 = self.recordheader.npw2
698
699
695
700 696 self.dataOut.channelList = range(1)
701 self.dataOut.nIncohInt = self.Num_inCoh
697 self.dataOut.nIncohInt = self.Num_inCoh
702 698 self.dataOut.nProfiles = self.Num_Bins
703 699 self.dataOut.nCohInt = 1
704 700 self.dataOut.windowOfFilter = 1
705 701 self.dataOut.utctime = dwell
706 self.dataOut.timeZone=0
707
702 self.dataOut.timeZone = 0
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
712
713 self.HSDVsign = numpy.fromfile( self.fp, [('HSDV',numpy.str_,4)],1)
714 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)
717
718 self.COFAsign = numpy.fromfile( self.fp, [('COFA',numpy.str_,4)],1)
719 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)
722
723 self.ZSPCsign = numpy.fromfile(self.fp, [('ZSPCsign',numpy.str_,4)],1)
724 self.SizeZSPC = numpy.fromfile(self.fp, [('SizeZSPC','<i4')],1)
725
726 self.dataOut.HSDV[0]=self.HSDV_Co[:][0]
727 self.dataOut.HSDV[1]=self.HSDV_Cx[:][0]
728
705 self.dataOut.heightList = self.SPARrawGate1 * self.__deltaHeigth + \
706 numpy.array(range(self.Num_Hei)) * self.__deltaHeigth
707
708 self.HSDVsign = numpy.fromfile(self.fp, [('HSDV', numpy.str_, 4)], 1)
709 self.SizeHSDV = numpy.fromfile(self.fp, [('SizeHSDV', '<i4')], 1)
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)
714
715 self.COFAsign = numpy.fromfile(self.fp, [('COFA', numpy.str_, 4)], 1)
716 self.SizeCOFA = numpy.fromfile(self.fp, [('SizeCOFA', '<i4')], 1)
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)
721
722 self.ZSPCsign = numpy.fromfile(
723 self.fp, [('ZSPCsign', numpy.str_, 4)], 1)
724 self.SizeZSPC = numpy.fromfile(self.fp, [('SizeZSPC', '<i4')], 1)
725
726 self.dataOut.HSDV[0] = self.HSDV_Co[:][0]
727 self.dataOut.HSDV[1] = self.HSDV_Cx[:][0]
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
731
730 # Number of spectral sub pieces containing significant power
731 nspc = numpy.fromfile(self.fp, [('nspc', 'int16')], 1)[0][0]
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
735
736 #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
740
741 self.data_spc[irg,binIndex:binIndex+nbins,0] = self.data_spc[irg,binIndex:binIndex+nbins,0]+jbin/65530.*jmax
742
743 #Cx_Channel
744 jbin = numpy.fromfile(self.fp, [('jbin','uint16')],nbins)[0][0]
745 jmax = numpy.fromfile(self.fp, [('jmax','float32')],1)[0][0]
746
747
748 self.data_spc[irg,binIndex:binIndex+nbins,1] = self.data_spc[irg,binIndex:binIndex+nbins,1]+jbin/65530.*jmax
749
750 for bin in range(self.Num_Bins):
751
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]
755
756
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
739
740 # Co_Channel
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
745
746 self.data_spc[irg, binIndex:binIndex + nbins, 0] = self.data_spc[irg,
747 binIndex:binIndex + nbins, 0] + jbin / 65530. * jmax
748
749 # Cx_Channel
750 jbin = numpy.fromfile(
751 self.fp, [('jbin', 'uint16')], nbins)[0][0]
752 jmax = numpy.fromfile(self.fp, [('jmax', 'float32')], 1)[0][0]
753
754 self.data_spc[irg, binIndex:binIndex + nbins, 1] = self.data_spc[irg,
755 binIndex:binIndex + nbins, 1] + jbin / 65530. * jmax
756
757 for bin in range(self.Num_Bins):
758
759 self.data_spc[:, bin, 0] = self.data_spc[:,
760 bin, 0] - self.dataOut.HSDV[:, 0]
761
762 self.data_spc[:, bin, 1] = self.data_spc[:,
763 bin, 1] - self.dataOut.HSDV[:, 1]
764
757 765 numpy.set_printoptions(threshold='nan')
758
759 self.data_spc = numpy.where(self.data_spc > 0. , self.data_spc, 0)
760
761 self.dataOut.COFA = numpy.array([self.COFA_Co , self.COFA_Cx])
762
766
767 self.data_spc = numpy.where(self.data_spc > 0., self.data_spc, 0)
768
769 self.dataOut.COFA = numpy.array([self.COFA_Co, self.COFA_Cx])
770
763 771 print ' '
764 print 'SPC',numpy.shape(self.dataOut.data_spc)
765 #print 'SPC',self.dataOut.data_spc
766
772 print 'SPC', numpy.shape(self.dataOut.data_spc)
773 # print 'SPC',self.dataOut.data_spc
774
767 775 noinor1 = 713031680
768 776 noinor2 = 30
769
770 npw1 = 1#0**(npw1/10) * noinor1 * noinor2
771 npw2 = 1#0**(npw2/10) * noinor1 * noinor2
777
778 npw1 = 1 # 0**(npw1/10) * noinor1 * noinor2
779 npw2 = 1 # 0**(npw2/10) * noinor1 * noinor2
772 780 self.dataOut.NPW = numpy.array([npw1, npw2])
773
781
774 782 print ' '
775
776 self.data_spc = numpy.transpose(self.data_spc, (2,1,0))
777 self.data_spc = numpy.fft.fftshift(self.data_spc, axes = 1)
778
783
784 self.data_spc = numpy.transpose(self.data_spc, (2, 1, 0))
785 self.data_spc = numpy.fft.fftshift(self.data_spc, axes=1)
786
779 787 self.data_spc = numpy.fliplr(self.data_spc)
780
781 self.data_spc = numpy.where(self.data_spc > 0. , self.data_spc, 0)
782 self.dataOut_spc= numpy.ones([1, self.Num_Bins , self.Num_Hei])
783 self.dataOut_spc[0,:,:] = self.data_spc[0,:,:]
784 #print 'SHAPE', self.dataOut_spc.shape
785 #For nyquist correction:
786 #fix = 20 # ~3m/s
788
789 self.data_spc = numpy.where(self.data_spc > 0., self.data_spc, 0)
790 self.dataOut_spc = numpy.ones([1, self.Num_Bins, self.Num_Hei])
791 self.dataOut_spc[0, :, :] = self.data_spc[0, :, :]
792 # print 'SHAPE', self.dataOut_spc.shape
793 # For nyquist correction:
794 # fix = 20 # ~3m/s
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
790
791
797
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
Show file before | Show file after | Hide comments
@@ -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