##// END OF EJS Templates
schain
RSS
Se modifico la forma de calcular la phase de la coherencia
Miguel Valdez -
r264:338ee2ca09fe
parent child
Show More
Show file before | Show file after | Hide comments
@@ -2,17 +2,16 import numpy
2 import datetime
2 import datetime
3 import sys
3 import sys
4 import matplotlib
4 import matplotlib
5 if sys.platform == 'linux':
5
6 matplotlib.use("GTKAgg")
6 if 'linux' in sys.platform:
7 if sys.platform == 'darwin':
7 matplotlib.use("TKAgg")
8
9 if 'darwin' in sys.platform:
8 matplotlib.use("TKAgg")
10 matplotlib.use("TKAgg")
11
9 import matplotlib.pyplot
12 import matplotlib.pyplot
10 import matplotlib.dates
11 #import scitools.numpyutils
12 from mpl_toolkits.axes_grid1 import make_axes_locatable
13
13
14 from matplotlib.dates import DayLocator, HourLocator, MinuteLocator, SecondLocator, DateFormatter
14 from mpl_toolkits.axes_grid1 import make_axes_locatable
15 from matplotlib.ticker import FuncFormatter
16 from matplotlib.ticker import *
15 from matplotlib.ticker import *
17
16
18 ###########################################
17 ###########################################
Show file before | Show file after | Hide comments
@@ -360,10 +360,10 class JRODataReader(JRODataIO, ProcessingUnit):
360 thisDate += datetime.timedelta(1)
360 thisDate += datetime.timedelta(1)
361
361
362 if pathList == []:
362 if pathList == []:
363 print "Any folder found into date range %s-%s" %(startDate, endDate)
363 print "Any folder was found for the date range: %s-%s" %(startDate, endDate)
364 return None, None
364 return None, None
365
365
366 print "%d folder(s) found [%s, ...]" %(len(pathList), pathList[0])
366 print "%d folder(s) was(were) found for the date range: %s-%s" %(len(pathList), startDate, endDate)
367
367
368 filenameList = []
368 filenameList = []
369 for thisPath in pathList:
369 for thisPath in pathList:
@@ -379,8 +379,10 class JRODataReader(JRODataIO, ProcessingUnit):
379 filenameList.append(filename)
379 filenameList.append(filename)
380
380
381 if not(filenameList):
381 if not(filenameList):
382 print "Any file found into time range %s-%s" %(startTime, endTime)
382 print "Any file was found for the time range %s - %s" %(startTime, endTime)
383 return None, None
383 return None, None
384
385 print "%d file(s) was(were) found for the time range: %s - %s" %(len(filenameList), startTime, endTime)
384
386
385 self.filenameList = filenameList
387 self.filenameList = filenameList
386
388
@@ -2077,10 +2079,19 class SpectraReader(JRODataReader):
2077
2079
2078 self.dataOut.flagShiftFFT = self.processingHeaderObj.shif_fft
2080 self.dataOut.flagShiftFFT = self.processingHeaderObj.shif_fft
2079
2081
2080 self.dataOut.flagDecodeData = True #asumo q la data no esta decodificada
2082 self.dataOut.flagDecodeData = False #asumo q la data no esta decodificada
2081
2083
2082 self.dataOut.flagDeflipData = True #asumo q la data no esta sin flip
2084 self.dataOut.flagDeflipData = True #asumo q la data no esta sin flip
2083
2085
2086 if self.processingHeaderObj.code != None:
2087
2088 self.dataOut.nCode = self.processingHeaderObj.nCode
2089
2090 self.dataOut.nBaud = self.processingHeaderObj.nBaud
2091
2092 self.dataOut.code = self.processingHeaderObj.code
2093
2094 self.dataOut.flagDecodeData = True
2084
2095
2085 return self.dataOut.data_spc
2096 return self.dataOut.data_spc
2086
2097
Show file before | Show file after | Hide comments
@@ -54,7 +54,8 class CrossSpectraPlot(Figure):
54
54
55 def run(self, dataOut, idfigure, wintitle="", pairsList=None, showprofile='True',
55 def run(self, dataOut, idfigure, wintitle="", pairsList=None, showprofile='True',
56 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
56 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
57 save=False, figpath='./', figfile=None):
57 save=False, figpath='./', figfile=None,
58 power_cmap='jet', coherence_cmap='jet', phase_cmap='RdBu_r'):
58
59
59 """
60 """
60
61
@@ -90,8 +91,8 class CrossSpectraPlot(Figure):
90 x = dataOut.getVelRange(1)
91 x = dataOut.getVelRange(1)
91 y = dataOut.getHeiRange()
92 y = dataOut.getHeiRange()
92 z = dataOut.data_spc[:,:,:]/factor
93 z = dataOut.data_spc[:,:,:]/factor
93 z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
94 # z = numpy.where(numpy.isfinite(z), z, numpy.NAN)
94 avg = numpy.average(numpy.abs(z), axis=1)
95 avg = numpy.abs(numpy.average(z, axis=1))
95 noise = dataOut.getNoise()/factor
96 noise = dataOut.getNoise()/factor
96
97
97 zdB = 10*numpy.log10(z)
98 zdB = 10*numpy.log10(z)
@@ -133,7 +134,7 class CrossSpectraPlot(Figure):
133 axes0.pcolor(x, y, zdB,
134 axes0.pcolor(x, y, zdB,
134 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
135 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
135 xlabel=xlabel, ylabel=ylabel, title=title,
136 xlabel=xlabel, ylabel=ylabel, title=title,
136 ticksize=9, cblabel='')
137 ticksize=9, colormap=power_cmap, cblabel='')
137
138
138 title = "Channel %d: %4.2fdB" %(pair[1], noisedB[pair[1]])
139 title = "Channel %d: %4.2fdB" %(pair[1], noisedB[pair[1]])
139 zdB = 10.*numpy.log10(dataOut.data_spc[pair[1],:,:]/factor)
140 zdB = 10.*numpy.log10(dataOut.data_spc[pair[1],:,:]/factor)
@@ -141,26 +142,26 class CrossSpectraPlot(Figure):
141 axes0.pcolor(x, y, zdB,
142 axes0.pcolor(x, y, zdB,
142 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
143 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax,
143 xlabel=xlabel, ylabel=ylabel, title=title,
144 xlabel=xlabel, ylabel=ylabel, title=title,
144 ticksize=9, cblabel='')
145 ticksize=9, colormap=power_cmap, cblabel='')
145
146
146 coherenceComplex = dataOut.data_cspc[pairsIndexList[i],:,:]/numpy.sqrt(dataOut.data_spc[pair[0],:,:]*dataOut.data_spc[pair[1],:,:])
147 coherenceComplex = dataOut.data_cspc[pairsIndexList[i],:,:]/numpy.sqrt(dataOut.data_spc[pair[0],:,:]*dataOut.data_spc[pair[1],:,:])
147 coherence = numpy.abs(coherenceComplex)
148 coherence = numpy.abs(coherenceComplex)
148 phase = numpy.arctan(-1*coherenceComplex.imag/coherenceComplex.real)*180/numpy.pi
149 # phase = numpy.arctan(-1*coherenceComplex.imag/coherenceComplex.real)*180/numpy.pi
149
150 phase = numpy.arctan2(coherenceComplex.imag, coherenceComplex.real)*180/numpy.pi
150
151
151 title = "Coherence %d%d" %(pair[0], pair[1])
152 title = "Coherence %d%d" %(pair[0], pair[1])
152 axes0 = self.axesList[i*self.__nsubplots+2]
153 axes0 = self.axesList[i*self.__nsubplots+2]
153 axes0.pcolor(x, y, coherence,
154 axes0.pcolor(x, y, coherence,
154 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=0, zmax=1,
155 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=0, zmax=1,
155 xlabel=xlabel, ylabel=ylabel, title=title,
156 xlabel=xlabel, ylabel=ylabel, title=title,
156 ticksize=9, cblabel='')
157 ticksize=9, colormap=coherence_cmap, cblabel='')
157
158
158 title = "Phase %d%d" %(pair[0], pair[1])
159 title = "Phase %d%d" %(pair[0], pair[1])
159 axes0 = self.axesList[i*self.__nsubplots+3]
160 axes0 = self.axesList[i*self.__nsubplots+3]
160 axes0.pcolor(x, y, phase,
161 axes0.pcolor(x, y, phase,
161 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=-180, zmax=180,
162 xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=-180, zmax=180,
162 xlabel=xlabel, ylabel=ylabel, title=title,
163 xlabel=xlabel, ylabel=ylabel, title=title,
163 ticksize=9, cblabel='', colormap='RdBu_r')
164 ticksize=9, colormap=phase_cmap, cblabel='')
164
165
165
166
166
167
@@ -719,7 +720,7 class CoherenceMap(Figure):
719
720
720 WIDTHPROF = None
721 WIDTHPROF = None
721 HEIGHTPROF = None
722 HEIGHTPROF = None
722 PREFIX = 'coherencemap'
723 PREFIX = 'cmap'
723
724
724 def __init__(self):
725 def __init__(self):
725 self.timerange = 2*60*60
726 self.timerange = 2*60*60
@@ -766,7 +767,8 class CoherenceMap(Figure):
766 def run(self, dataOut, idfigure, wintitle="", pairsList=None, showprofile='True',
767 def run(self, dataOut, idfigure, wintitle="", pairsList=None, showprofile='True',
767 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
768 xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None,
768 timerange=None,
769 timerange=None,
769 save=False, figpath='./', figfile=None):
770 save=False, figpath='./', figfile=None,
771 coherence_cmap='jet', phase_cmap='RdBu_r'):
770
772
771 if pairsList == None:
773 if pairsList == None:
772 pairsIndexList = dataOut.pairsIndexList
774 pairsIndexList = dataOut.pairsIndexList
@@ -817,9 +819,12 class CoherenceMap(Figure):
817
819
818 pair = dataOut.pairsList[pairsIndexList[i]]
820 pair = dataOut.pairsList[pairsIndexList[i]]
819 coherenceComplex = dataOut.data_cspc[pairsIndexList[i],:,:]/numpy.sqrt(dataOut.data_spc[pair[0],:,:]*dataOut.data_spc[pair[1],:,:])
821 coherenceComplex = dataOut.data_cspc[pairsIndexList[i],:,:]/numpy.sqrt(dataOut.data_spc[pair[0],:,:]*dataOut.data_spc[pair[1],:,:])
820 coherence = numpy.abs(coherenceComplex)
822 avgcoherenceComplex = numpy.average(coherenceComplex, axis=0)
821 avg = numpy.average(coherence, axis=0)
823 coherence = numpy.abs(avgcoherenceComplex)
822 z = avg.reshape((1,-1))
824 # coherence = numpy.abs(coherenceComplex)
825 # avg = numpy.average(coherence, axis=0)
826
827 z = coherence.reshape((1,-1))
823
828
824 counter = 0
829 counter = 0
825
830
@@ -828,33 +833,34 class CoherenceMap(Figure):
828 axes.pcolor(x, y, z,
833 axes.pcolor(x, y, z,
829 xmin=tmin, xmax=tmax, ymin=ymin, ymax=ymax, zmin=0, zmax=1,
834 xmin=tmin, xmax=tmax, ymin=ymin, ymax=ymax, zmin=0, zmax=1,
830 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
835 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
831 ticksize=9, cblabel='', cbsize="1%")
836 ticksize=9, cblabel='', colormap=coherence_cmap, cbsize="1%")
832
837
833 if self.__showprofile:
838 if self.__showprofile:
834 counter += 1
839 counter += 1
835 axes = self.axesList[i*self.__nsubplots*2 + counter]
840 axes = self.axesList[i*self.__nsubplots*2 + counter]
836 axes.pline(avg, y,
841 axes.pline(coherence, y,
837 xmin=0, xmax=1, ymin=ymin, ymax=ymax,
842 xmin=0, xmax=1, ymin=ymin, ymax=ymax,
838 xlabel='', ylabel='', title='', ticksize=7,
843 xlabel='', ylabel='', title='', ticksize=7,
839 ytick_visible=False, nxticks=5,
844 ytick_visible=False, nxticks=5,
840 grid='x')
845 grid='x')
841
846
842 counter += 1
847 counter += 1
843 phase = numpy.arctan(-1*coherenceComplex.imag/coherenceComplex.real)*180/numpy.pi
848 # phase = numpy.arctan(-1*coherenceComplex.imag/coherenceComplex.real)*180/numpy.pi
844 avg = numpy.average(phase, axis=0)
849 phase = numpy.arctan2(avgcoherenceComplex.imag, avgcoherenceComplex.real)*180/numpy.pi
845 z = avg.reshape((1,-1))
850 # avg = numpy.average(phase, axis=0)
851 z = phase.reshape((1,-1))
846
852
847 title = "Phase %d%d: %s" %(pair[0], pair[1], thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
853 title = "Phase %d%d: %s" %(pair[0], pair[1], thisDatetime.strftime("%d-%b-%Y %H:%M:%S"))
848 axes = self.axesList[i*self.__nsubplots*2 + counter]
854 axes = self.axesList[i*self.__nsubplots*2 + counter]
849 axes.pcolor(x, y, z,
855 axes.pcolor(x, y, z,
850 xmin=tmin, xmax=tmax, ymin=ymin, ymax=ymax, zmin=-180, zmax=180,
856 xmin=tmin, xmax=tmax, ymin=ymin, ymax=ymax, zmin=-180, zmax=180,
851 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
857 xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True,
852 ticksize=9, cblabel='', colormap='RdBu', cbsize="1%")
858 ticksize=9, cblabel='', colormap=phase_cmap, cbsize="1%")
853
859
854 if self.__showprofile:
860 if self.__showprofile:
855 counter += 1
861 counter += 1
856 axes = self.axesList[i*self.__nsubplots*2 + counter]
862 axes = self.axesList[i*self.__nsubplots*2 + counter]
857 axes.pline(avg, y,
863 axes.pline(phase, y,
858 xmin=-180, xmax=180, ymin=ymin, ymax=ymax,
864 xmin=-180, xmax=180, ymin=ymin, ymax=ymax,
859 xlabel='', ylabel='', title='', ticksize=7,
865 xlabel='', ylabel='', title='', ticksize=7,
860 ytick_visible=False, nxticks=4,
866 ytick_visible=False, nxticks=4,
General Comments 0
You need to be logged in to leave comments. Login now