@@ -0,0 +1,135 | |||
|
1 | # DIAS 19 Y 20 FEB 2014 | |
|
2 | # Comprobacion de Resultados DBS con SA | |
|
3 | ||
|
4 | import os, sys | |
|
5 | ||
|
6 | path = os.path.split(os.getcwd())[0] | |
|
7 | sys.path.append(path) | |
|
8 | ||
|
9 | from controller import * | |
|
10 | ||
|
11 | desc = "DBS Experiment Test" | |
|
12 | filename = "DBStest.xml" | |
|
13 | ||
|
14 | controllerObj = Project() | |
|
15 | ||
|
16 | controllerObj.setup(id = '191', name='test01', description=desc) | |
|
17 | ||
|
18 | #Experimentos | |
|
19 | ||
|
20 | path = '/host/Jicamarca/EW_Drifts/d2012248' | |
|
21 | pathFigure = '/home/propietario/workspace/Graficos/drifts' | |
|
22 | ||
|
23 | ||
|
24 | path = "/home/soporte/Data/drifts" | |
|
25 | pathFigure = '/home/soporte/workspace/Graficos/drifts/prueba' | |
|
26 | ||
|
27 | xmin = 11.75 | |
|
28 | xmax = 14.75 | |
|
29 | #------------------------------------------------------------------------------------------------ | |
|
30 | readUnitConfObj = controllerObj.addReadUnit(datatype='VoltageReader', | |
|
31 | path=path, | |
|
32 | startDate='2012/01/01', | |
|
33 | endDate='2012/12/31', | |
|
34 | startTime='00:00:00', | |
|
35 | endTime='23:59:59', | |
|
36 | online=0, | |
|
37 | walk=1) | |
|
38 | ||
|
39 | opObj11 = readUnitConfObj.addOperation(name='printNumberOfBlock') | |
|
40 | ||
|
41 | #-------------------------------------------------------------------------------------------------- | |
|
42 | ||
|
43 | procUnitConfObj0 = controllerObj.addProcUnit(datatype='VoltageProc', inputId=readUnitConfObj.getId()) | |
|
44 | ||
|
45 | opObj11 = procUnitConfObj0.addOperation(name='ProfileSelector', optype='other') | |
|
46 | opObj11.addParameter(name='profileRangeList', value='0,127', format='intlist') | |
|
47 | ||
|
48 | opObj11 = procUnitConfObj0.addOperation(name='filterByHeights') | |
|
49 | opObj11.addParameter(name='window', value='3', format='int') | |
|
50 | ||
|
51 | opObj11 = procUnitConfObj0.addOperation(name='Decoder', optype='other') | |
|
52 | # opObj11.addParameter(name='code', value='1,-1', format='floatlist') | |
|
53 | # opObj11.addParameter(name='nCode', value='2', format='int') | |
|
54 | # opObj11.addParameter(name='nBaud', value='1', format='int') | |
|
55 | ||
|
56 | procUnitConfObj1 = controllerObj.addProcUnit(datatype='SpectraProc', inputId=procUnitConfObj0.getId()) | |
|
57 | procUnitConfObj1.addParameter(name='nFFTPoints', value='128', format='int') | |
|
58 | procUnitConfObj1.addParameter(name='nProfiles', value='128', format='int') | |
|
59 | procUnitConfObj1.addParameter(name='pairsList', value='(0,1),(2,3)', format='pairsList')#,(2,3) | |
|
60 | ||
|
61 | opObj11 = procUnitConfObj1.addOperation(name='selectHeights') | |
|
62 | # # opObj11.addParameter(name='minHei', value='320.0', format='float') | |
|
63 | # # opObj11.addParameter(name='maxHei', value='350.0', format='float') | |
|
64 | opObj11.addParameter(name='minHei', value='200.0', format='float') | |
|
65 | opObj11.addParameter(name='maxHei', value='600.0', format='float') | |
|
66 | ||
|
67 | opObj11 = procUnitConfObj1.addOperation(name='selectChannels') | |
|
68 | opObj11.addParameter(name='channelList', value='0,1,2,3', format='intlist') | |
|
69 | ||
|
70 | opObj11 = procUnitConfObj1.addOperation(name='IncohInt', optype='other') | |
|
71 | opObj11.addParameter(name='timeInterval', value='300.0', format='float') | |
|
72 | ||
|
73 | opObj13 = procUnitConfObj1.addOperation(name='removeDC') | |
|
74 | ||
|
75 | # opObj14 = procUnitConfObj1.addOperation(name='SpectraPlot', optype='other') | |
|
76 | # opObj14.addParameter(name='id', value='1', format='int') | |
|
77 | # # opObj14.addParameter(name='wintitle', value='Con interf', format='str') | |
|
78 | # opObj14.addParameter(name='save', value='1', format='bool') | |
|
79 | # opObj14.addParameter(name='figpath', value=pathFigure, format='str') | |
|
80 | # # opObj14.addParameter(name='zmin', value='5', format='int') | |
|
81 | # opObj14.addParameter(name='zmax', value='30', format='int') | |
|
82 | # | |
|
83 | # opObj12 = procUnitConfObj1.addOperation(name='RTIPlot', optype='other') | |
|
84 | # opObj12.addParameter(name='id', value='2', format='int') | |
|
85 | # opObj12.addParameter(name='wintitle', value='RTI Plot', format='str') | |
|
86 | # opObj12.addParameter(name='save', value='1', format='bool') | |
|
87 | # opObj12.addParameter(name='figpath', value = pathFigure, format='str') | |
|
88 | # opObj12.addParameter(name='xmin', value=xmin, format='float') | |
|
89 | # opObj12.addParameter(name='xmax', value=xmax, format='float') | |
|
90 | # # opObj12.addParameter(name='zmin', value='5', format='int') | |
|
91 | # opObj12.addParameter(name='zmax', value='30', format='int') | |
|
92 | ||
|
93 | #-------------------------------------------------------------------------------------------------- | |
|
94 | ||
|
95 | procUnitConfObj2 = controllerObj.addProcUnit(datatype='ParametersProc', inputId=procUnitConfObj1.getId()) | |
|
96 | opObj20 = procUnitConfObj2.addOperation(name='SpectralFitting') | |
|
97 | opObj20.addParameter(name='path', value='/home/soporte/workspace/RemoteSystemsTempFiles', format='str') | |
|
98 | opObj20.addParameter(name='file', value='modelSpectralFitting', format='str') | |
|
99 | opObj20.addParameter(name='groupList', value='(0,1),(2,3)',format='multiList') | |
|
100 | ||
|
101 | opObj11 = procUnitConfObj2.addOperation(name='SpectralFittingPlot', optype='other') | |
|
102 | opObj11.addParameter(name='id', value='3', format='int') | |
|
103 | opObj11.addParameter(name='wintitle', value='DopplerPlot', format='str') | |
|
104 | opObj11.addParameter(name='cutHeight', value='350', format='int') | |
|
105 | opObj11.addParameter(name='fit', value='1', format='int')#1--True/include fit | |
|
106 | opObj11.addParameter(name='save', value='1', format='bool') | |
|
107 | opObj11.addParameter(name='figpath', value = pathFigure, format='str') | |
|
108 | ||
|
109 | opObj11 = procUnitConfObj2.addOperation(name='EWDriftsEstimation', optype='other') | |
|
110 | opObj11.addParameter(name='zenith', value='-3.80208,3.10658', format='floatlist') | |
|
111 | opObj11.addParameter(name='zenithCorrection', value='0.183201', format='float') | |
|
112 | ||
|
113 | opObj23 = procUnitConfObj2.addOperation(name='EWDriftsPlot', optype='other') | |
|
114 | opObj23.addParameter(name='id', value='4', format='int') | |
|
115 | opObj23.addParameter(name='wintitle', value='EW Drifts', format='str') | |
|
116 | opObj23.addParameter(name='save', value='1', format='bool') | |
|
117 | opObj23.addParameter(name='figpath', value = pathFigure, format='str') | |
|
118 | opObj23.addParameter(name='zminZonal', value='-150', format='int') | |
|
119 | opObj23.addParameter(name='zmaxZonal', value='150', format='int') | |
|
120 | opObj23.addParameter(name='zminVertical', value='-30', format='float') | |
|
121 | opObj23.addParameter(name='zmaxVertical', value='30', format='float') | |
|
122 | opObj23.addParameter(name='SNR_1', value='1', format='bool') | |
|
123 | opObj23.addParameter(name='SNRmax', value='5', format='int') | |
|
124 | # opObj23.addParameter(name='SNRthresh', value='-50', format='float') | |
|
125 | opObj23.addParameter(name='xmin', value=xmin, format='float') | |
|
126 | opObj23.addParameter(name='xmax', value=xmax, format='float') | |
|
127 | #-------------------------------------------------------------------------------------------------- | |
|
128 | print "Escribiendo el archivo XML" | |
|
129 | controllerObj.writeXml(filename) | |
|
130 | print "Leyendo el archivo XML" | |
|
131 | controllerObj.readXml(filename) | |
|
132 | ||
|
133 | controllerObj.createObjects() | |
|
134 | controllerObj.connectObjects() | |
|
135 | controllerObj.run() No newline at end of file |
@@ -917,10 +917,18 class Correlation(JROData): | |||
|
917 | 917 | |
|
918 | 918 | class Parameters(JROData): |
|
919 | 919 | |
|
920 | #Information from previous data | |
|
921 | ||
|
920 | 922 | inputUnit = None #Type of data to be processed |
|
921 | 923 | |
|
922 | 924 | operation = None #Type of operation to parametrize |
|
923 | 925 | |
|
926 | normFactor = None #Normalization Factor | |
|
927 | ||
|
928 | groupList = None #List of Pairs, Groups, etc | |
|
929 | ||
|
930 | #Parameters | |
|
931 | ||
|
924 | 932 | data_param = None #Parameters obtained |
|
925 | 933 | |
|
926 | 934 | data_pre = None #Data Pre Parametrization |
@@ -933,17 +941,25 class Parameters(JROData): | |||
|
933 | 941 | |
|
934 | 942 | SNR = None #Signal to Noise Ratio |
|
935 | 943 | |
|
936 | pairsList = None #List of Pairs for Cross correlations or Cross spectrum | |
|
937 | ||
|
938 | 944 | initUtcTime = None #Initial UTC time |
|
939 | 945 | |
|
940 | 946 | paramInterval = None #Time interval to calculate Parameters in seconds |
|
941 | 947 | |
|
942 | windsInterval = None #Time interval to calculate Winds in seconds | |
|
948 | #Fitting | |
|
949 | ||
|
950 | constants = None | |
|
951 | ||
|
952 | error = None | |
|
953 | ||
|
954 | library = None | |
|
955 | ||
|
956 | #Output signal | |
|
957 | ||
|
958 | outputInterval = None #Time interval to calculate output signal in seconds | |
|
959 | ||
|
960 | data_output = None #Out signal | |
|
943 | 961 | |
|
944 | normFactor = None #Normalization Factor | |
|
945 | 962 | |
|
946 | winds = None #Wind estimations | |
|
947 | 963 | |
|
948 | 964 | def __init__(self): |
|
949 | 965 | ''' |
@@ -960,7 +976,7 class Parameters(JROData): | |||
|
960 | 976 | datatime = [] |
|
961 | 977 | |
|
962 | 978 | datatime.append(self.initUtcTime) |
|
963 |
datatime.append(self.initUtcTime + self. |
|
|
979 | datatime.append(self.initUtcTime + self.outputInterval - 1) | |
|
964 | 980 | |
|
965 | 981 | datatime = numpy.array(datatime) |
|
966 | 982 |
@@ -455,7 +455,7 class WindProfilerPlot(Figure): | |||
|
455 | 455 | # y = dataOut.heightRange |
|
456 | 456 | y = dataOut.heightRange |
|
457 | 457 | |
|
458 |
z = dataOut. |
|
|
458 | z = dataOut.data_output.copy() | |
|
459 | 459 | nplots = z.shape[0] #Number of wind dimensions estimated |
|
460 | 460 | nplotsw = nplots |
|
461 | 461 | |
@@ -771,4 +771,408 class ParametersPlot(Figure): | |||
|
771 | 771 | if x[1] >= self.axesList[0].xmax: |
|
772 | 772 | self.counter_imagwr = wr_period |
|
773 | 773 | self.__isConfig = False |
|
774 | self.figfile = None | |
|
775 | ||
|
776 | ||
|
777 | class SpectralFittingPlot(Figure): | |
|
778 | ||
|
779 | __isConfig = None | |
|
780 | __nsubplots = None | |
|
781 | ||
|
782 | WIDTHPROF = None | |
|
783 | HEIGHTPROF = None | |
|
784 | PREFIX = 'prm' | |
|
785 | ||
|
786 | ||
|
787 | N = None | |
|
788 | ippSeconds = None | |
|
789 | ||
|
790 | def __init__(self): | |
|
791 | self.__isConfig = False | |
|
792 | self.__nsubplots = 1 | |
|
793 | ||
|
794 | self.WIDTH = 450 | |
|
795 | self.HEIGHT = 250 | |
|
796 | self.WIDTHPROF = 0 | |
|
797 | self.HEIGHTPROF = 0 | |
|
798 | ||
|
799 | def getSubplots(self): | |
|
800 | ||
|
801 | ncol = int(numpy.sqrt(self.nplots)+0.9) | |
|
802 | nrow = int(self.nplots*1./ncol + 0.9) | |
|
803 | ||
|
804 | return nrow, ncol | |
|
805 | ||
|
806 | def setup(self, id, nplots, wintitle, showprofile=False, show=True): | |
|
807 | ||
|
808 | showprofile = False | |
|
809 | self.__showprofile = showprofile | |
|
810 | self.nplots = nplots | |
|
811 | ||
|
812 | ncolspan = 5 | |
|
813 | colspan = 4 | |
|
814 | if showprofile: | |
|
815 | ncolspan = 5 | |
|
816 | colspan = 4 | |
|
817 | self.__nsubplots = 2 | |
|
818 | ||
|
819 | self.createFigure(id = id, | |
|
820 | wintitle = wintitle, | |
|
821 | widthplot = self.WIDTH + self.WIDTHPROF, | |
|
822 | heightplot = self.HEIGHT + self.HEIGHTPROF, | |
|
823 | show=show) | |
|
824 | ||
|
825 | nrow, ncol = self.getSubplots() | |
|
826 | ||
|
827 | counter = 0 | |
|
828 | for y in range(nrow): | |
|
829 | for x in range(ncol): | |
|
830 | ||
|
831 | if counter >= self.nplots: | |
|
832 | break | |
|
833 | ||
|
834 | self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan, colspan, 1) | |
|
835 | ||
|
836 | if showprofile: | |
|
837 | self.addAxes(nrow, ncol*ncolspan, y, x*ncolspan+colspan, 1, 1) | |
|
838 | ||
|
839 | counter += 1 | |
|
840 | ||
|
841 | def run(self, dataOut, id, cutHeight=None, fit=False, wintitle="", channelList=None, showprofile=True, | |
|
842 | xmin=None, xmax=None, ymin=None, ymax=None, | |
|
843 | save=False, figpath='./', figfile=None, show=True): | |
|
844 | ||
|
845 | """ | |
|
846 | ||
|
847 | Input: | |
|
848 | dataOut : | |
|
849 | id : | |
|
850 | wintitle : | |
|
851 | channelList : | |
|
852 | showProfile : | |
|
853 | xmin : None, | |
|
854 | xmax : None, | |
|
855 | zmin : None, | |
|
856 | zmax : None | |
|
857 | """ | |
|
858 | ||
|
859 | if cutHeight==None: | |
|
860 | h=270 | |
|
861 | heightindex = numpy.abs(cutHeight - dataOut.heightList).argmin() | |
|
862 | cutHeight = dataOut.heightList[heightindex] | |
|
863 | ||
|
864 | factor = dataOut.normFactor | |
|
865 | x = dataOut.abscissaRange[:-1] | |
|
866 | #y = dataOut.getHeiRange() | |
|
867 | ||
|
868 | z = dataOut.data_pre[:,:,heightindex]/factor | |
|
869 | z = numpy.where(numpy.isfinite(z), z, numpy.NAN) | |
|
870 | avg = numpy.average(z, axis=1) | |
|
871 | listChannels = z.shape[0] | |
|
872 | ||
|
873 | #Reconstruct Function | |
|
874 | if fit==True: | |
|
875 | groupArray = dataOut.groupList | |
|
876 | listChannels = groupArray.reshape((groupArray.size)) | |
|
877 | listChannels.sort() | |
|
878 | spcFitLine = numpy.zeros(z.shape) | |
|
879 | constants = dataOut.constants | |
|
880 | ||
|
881 | nGroups = groupArray.shape[0] | |
|
882 | nChannels = groupArray.shape[1] | |
|
883 | nProfiles = z.shape[1] | |
|
884 | ||
|
885 | for f in range(nGroups): | |
|
886 | groupChann = groupArray[f,:] | |
|
887 | p = dataOut.data_param[f,:,heightindex] | |
|
888 | # p = numpy.array([ 89.343967,0.14036615,0.17086219,18.89835291,1.58388365,1.55099167]) | |
|
889 | fitLineAux = dataOut.library.modelFunction(p, constants)*nProfiles | |
|
890 | fitLineAux = fitLineAux.reshape((nChannels,nProfiles)) | |
|
891 | spcFitLine[groupChann,:] = fitLineAux | |
|
892 | # spcFitLine = spcFitLine/factor | |
|
893 | ||
|
894 | z = z[listChannels,:] | |
|
895 | spcFitLine = spcFitLine[listChannels,:] | |
|
896 | spcFitLinedB = 10*numpy.log10(spcFitLine) | |
|
897 | ||
|
898 | zdB = 10*numpy.log10(z) | |
|
899 | #thisDatetime = dataOut.datatime | |
|
900 | thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[1]) | |
|
901 | title = wintitle + " Doppler Spectra: %s" %(thisDatetime.strftime("%d-%b-%Y %H:%M:%S")) | |
|
902 | xlabel = "Velocity (m/s)" | |
|
903 | ylabel = "Spectrum" | |
|
904 | ||
|
905 | if not self.__isConfig: | |
|
906 | ||
|
907 | nplots = listChannels.size | |
|
908 | ||
|
909 | self.setup(id=id, | |
|
910 | nplots=nplots, | |
|
911 | wintitle=wintitle, | |
|
912 | showprofile=showprofile, | |
|
913 | show=show) | |
|
914 | ||
|
915 | if xmin == None: xmin = numpy.nanmin(x) | |
|
916 | if xmax == None: xmax = numpy.nanmax(x) | |
|
917 | if ymin == None: ymin = numpy.nanmin(zdB) | |
|
918 | if ymax == None: ymax = numpy.nanmax(zdB)+2 | |
|
919 | ||
|
920 | self.__isConfig = True | |
|
921 | ||
|
922 | self.setWinTitle(title) | |
|
923 | for i in range(self.nplots): | |
|
924 | # title = "Channel %d: %4.2fdB" %(dataOut.channelList[i]+1, noisedB[i]) | |
|
925 | title = "Height %4.1f km\nChannel %d:" %(cutHeight, listChannels[i]+1) | |
|
926 | axes = self.axesList[i*self.__nsubplots] | |
|
927 | if fit == False: | |
|
928 | axes.pline(x, zdB[i,:], | |
|
929 | xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, | |
|
930 | xlabel=xlabel, ylabel=ylabel, title=title | |
|
931 | ) | |
|
932 | if fit == True: | |
|
933 | fitline=spcFitLinedB[i,:] | |
|
934 | y=numpy.vstack([zdB[i,:],fitline] ) | |
|
935 | legendlabels=['Data','Fitting'] | |
|
936 | axes.pmultilineyaxis(x, y, | |
|
937 | xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, | |
|
938 | xlabel=xlabel, ylabel=ylabel, title=title, | |
|
939 | legendlabels=legendlabels, marker=None, | |
|
940 | linestyle='solid', grid='both') | |
|
941 | ||
|
942 | self.draw() | |
|
943 | ||
|
944 | if save: | |
|
945 | date = thisDatetime.strftime("%Y%m%d_%H%M%S") | |
|
946 | if figfile == None: | |
|
947 | figfile = self.getFilename(name = date) | |
|
948 | ||
|
949 | self.saveFigure(figpath, figfile) | |
|
950 | ||
|
951 | ||
|
952 | class EWDriftsPlot(Figure): | |
|
953 | ||
|
954 | __isConfig = None | |
|
955 | __nsubplots = None | |
|
956 | ||
|
957 | WIDTHPROF = None | |
|
958 | HEIGHTPROF = None | |
|
959 | PREFIX = 'drift' | |
|
960 | ||
|
961 | def __init__(self): | |
|
962 | ||
|
963 | self.timerange = 2*60*60 | |
|
964 | self.isConfig = False | |
|
965 | self.__nsubplots = 1 | |
|
966 | ||
|
967 | self.WIDTH = 800 | |
|
968 | self.HEIGHT = 150 | |
|
969 | self.WIDTHPROF = 120 | |
|
970 | self.HEIGHTPROF = 0 | |
|
971 | self.counter_imagwr = 0 | |
|
972 | ||
|
973 | self.PLOT_CODE = 0 | |
|
974 | self.FTP_WEI = None | |
|
975 | self.EXP_CODE = None | |
|
976 | self.SUB_EXP_CODE = None | |
|
977 | self.PLOT_POS = None | |
|
978 | self.tmin = None | |
|
979 | self.tmax = None | |
|
980 | ||
|
981 | self.xmin = None | |
|
982 | self.xmax = None | |
|
983 | ||
|
984 | self.figfile = None | |
|
985 | ||
|
986 | def getSubplots(self): | |
|
987 | ||
|
988 | ncol = 1 | |
|
989 | nrow = self.nplots | |
|
990 | ||
|
991 | return nrow, ncol | |
|
992 | ||
|
993 | def setup(self, id, nplots, wintitle, showprofile=True, show=True): | |
|
994 | ||
|
995 | self.__showprofile = showprofile | |
|
996 | self.nplots = nplots | |
|
997 | ||
|
998 | ncolspan = 1 | |
|
999 | colspan = 1 | |
|
1000 | ||
|
1001 | self.createFigure(id = id, | |
|
1002 | wintitle = wintitle, | |
|
1003 | widthplot = self.WIDTH + self.WIDTHPROF, | |
|
1004 | heightplot = self.HEIGHT + self.HEIGHTPROF, | |
|
1005 | show=show) | |
|
1006 | ||
|
1007 | nrow, ncol = self.getSubplots() | |
|
1008 | ||
|
1009 | counter = 0 | |
|
1010 | for y in range(nrow): | |
|
1011 | if counter >= self.nplots: | |
|
1012 | break | |
|
1013 | ||
|
1014 | self.addAxes(nrow, ncol*ncolspan, y, 0, colspan, 1) | |
|
1015 | counter += 1 | |
|
1016 | ||
|
1017 | def run(self, dataOut, id, wintitle="", channelList=None, | |
|
1018 | xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None, | |
|
1019 | zmaxVertical = None, zminVertical = None, zmaxZonal = None, zminZonal = None, | |
|
1020 | timerange=None, SNRthresh = -numpy.inf, SNRmin = None, SNRmax = None, SNR_1 = False, | |
|
1021 | save=False, figpath='', lastone=0,figfile=None, ftp=False, wr_period=1, show=True, | |
|
1022 | server=None, folder=None, username=None, password=None, | |
|
1023 | ftp_wei=0, exp_code=0, sub_exp_code=0, plot_pos=0): | |
|
1024 | """ | |
|
1025 | ||
|
1026 | Input: | |
|
1027 | dataOut : | |
|
1028 | id : | |
|
1029 | wintitle : | |
|
1030 | channelList : | |
|
1031 | showProfile : | |
|
1032 | xmin : None, | |
|
1033 | xmax : None, | |
|
1034 | ymin : None, | |
|
1035 | ymax : None, | |
|
1036 | zmin : None, | |
|
1037 | zmax : None | |
|
1038 | """ | |
|
1039 | ||
|
1040 | if channelList == None: | |
|
1041 | channelIndexList = dataOut.channelIndexList | |
|
1042 | else: | |
|
1043 | channelIndexList = [] | |
|
1044 | for channel in channelList: | |
|
1045 | if channel not in dataOut.channelList: | |
|
1046 | raise ValueError, "Channel %d is not in dataOut.channelList" | |
|
1047 | channelIndexList.append(dataOut.channelList.index(channel)) | |
|
1048 | ||
|
1049 | if timerange != None: | |
|
1050 | self.timerange = timerange | |
|
1051 | ||
|
1052 | tmin = None | |
|
1053 | tmax = None | |
|
1054 | ||
|
1055 | x = dataOut.getTimeRange1() | |
|
1056 | # y = dataOut.heightRange | |
|
1057 | y = dataOut.heightList | |
|
1058 | ||
|
1059 | z = dataOut.data_output | |
|
1060 | nplots = z.shape[0] #Number of wind dimensions estimated | |
|
1061 | nplotsw = nplots | |
|
1062 | ||
|
1063 | #If there is a SNR function defined | |
|
1064 | if dataOut.SNR != None: | |
|
1065 | nplots += 1 | |
|
1066 | SNR = dataOut.SNR | |
|
1067 | ||
|
1068 | if SNR_1: | |
|
1069 | SNR += 1 | |
|
1070 | ||
|
1071 | SNRavg = numpy.average(SNR, axis=0) | |
|
1072 | ||
|
1073 | SNRdB = 10*numpy.log10(SNR) | |
|
1074 | SNRavgdB = 10*numpy.log10(SNRavg) | |
|
1075 | ||
|
1076 | ind = numpy.where(SNRavg < 10**(SNRthresh/10))[0] | |
|
1077 | ||
|
1078 | for i in range(nplotsw): | |
|
1079 | z[i,ind] = numpy.nan | |
|
1080 | ||
|
1081 | ||
|
1082 | showprofile = False | |
|
1083 | # thisDatetime = dataOut.datatime | |
|
1084 | thisDatetime = datetime.datetime.utcfromtimestamp(dataOut.getTimeRange()[1]) | |
|
1085 | title = wintitle + " EW Drifts" | |
|
1086 | xlabel = "" | |
|
1087 | ylabel = "Height (Km)" | |
|
1088 | ||
|
1089 | if not self.__isConfig: | |
|
1090 | ||
|
1091 | self.setup(id=id, | |
|
1092 | nplots=nplots, | |
|
1093 | wintitle=wintitle, | |
|
1094 | showprofile=showprofile, | |
|
1095 | show=show) | |
|
1096 | ||
|
1097 | self.xmin, self.xmax = self.getTimeLim(x, xmin, xmax, timerange) | |
|
1098 | ||
|
1099 | if ymin == None: ymin = numpy.nanmin(y) | |
|
1100 | if ymax == None: ymax = numpy.nanmax(y) | |
|
1101 | ||
|
1102 | if zmaxZonal == None: zmaxZonal = numpy.nanmax(abs(z[0,:])) | |
|
1103 | if zminZonal == None: zminZonal = -zmaxZonal | |
|
1104 | if zmaxVertical == None: zmaxVertical = numpy.nanmax(abs(z[1,:])) | |
|
1105 | if zminVertical == None: zminVertical = -zmaxVertical | |
|
1106 | ||
|
1107 | if dataOut.SNR != None: | |
|
1108 | if SNRmin == None: SNRmin = numpy.nanmin(SNRavgdB) | |
|
1109 | if SNRmax == None: SNRmax = numpy.nanmax(SNRavgdB) | |
|
1110 | ||
|
1111 | self.FTP_WEI = ftp_wei | |
|
1112 | self.EXP_CODE = exp_code | |
|
1113 | self.SUB_EXP_CODE = sub_exp_code | |
|
1114 | self.PLOT_POS = plot_pos | |
|
1115 | ||
|
1116 | self.name = thisDatetime.strftime("%Y%m%d_%H%M%S") | |
|
1117 | self.__isConfig = True | |
|
1118 | ||
|
1119 | ||
|
1120 | self.setWinTitle(title) | |
|
1121 | ||
|
1122 | if ((self.xmax - x[1]) < (x[1]-x[0])): | |
|
1123 | x[1] = self.xmax | |
|
1124 | ||
|
1125 | strWind = ['Zonal','Vertical'] | |
|
1126 | strCb = 'Velocity (m/s)' | |
|
1127 | zmaxVector = [zmaxZonal, zmaxVertical] | |
|
1128 | zminVector = [zminZonal, zminVertical] | |
|
1129 | ||
|
1130 | for i in range(nplotsw): | |
|
1131 | ||
|
1132 | title = "%s Drifts: %s" %(strWind[i], thisDatetime.strftime("%Y/%m/%d %H:%M:%S")) | |
|
1133 | axes = self.axesList[i*self.__nsubplots] | |
|
1134 | ||
|
1135 | z1 = z[i,:].reshape((1,-1)) | |
|
1136 | ||
|
1137 | axes.pcolorbuffer(x, y, z1, | |
|
1138 | xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=zminVector[i], zmax=zmaxVector[i], | |
|
1139 | xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True, | |
|
1140 | ticksize=9, cblabel=strCb, cbsize="1%", colormap="RdBu_r") | |
|
1141 | ||
|
1142 | if dataOut.SNR != None: | |
|
1143 | i += 1 | |
|
1144 | if SNR_1: | |
|
1145 | title = "Signal Noise Ratio + 1 (SNR+1): %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S")) | |
|
1146 | else: | |
|
1147 | title = "Signal Noise Ratio (SNR): %s" %(thisDatetime.strftime("%Y/%m/%d %H:%M:%S")) | |
|
1148 | axes = self.axesList[i*self.__nsubplots] | |
|
1149 | SNRavgdB = SNRavgdB.reshape((1,-1)) | |
|
1150 | ||
|
1151 | axes.pcolorbuffer(x, y, SNRavgdB, | |
|
1152 | xmin=self.xmin, xmax=self.xmax, ymin=ymin, ymax=ymax, zmin=SNRmin, zmax=SNRmax, | |
|
1153 | xlabel=xlabel, ylabel=ylabel, title=title, rti=True, XAxisAsTime=True, | |
|
1154 | ticksize=9, cblabel='', cbsize="1%", colormap="jet") | |
|
1155 | ||
|
1156 | self.draw() | |
|
1157 | ||
|
1158 | if self.figfile == None: | |
|
1159 | str_datetime = thisDatetime.strftime("%Y%m%d_%H%M%S") | |
|
1160 | self.figfile = self.getFilename(name = str_datetime) | |
|
1161 | ||
|
1162 | if figpath != '': | |
|
1163 | ||
|
1164 | self.counter_imagwr += 1 | |
|
1165 | if (self.counter_imagwr>=wr_period): | |
|
1166 | # store png plot to local folder | |
|
1167 | self.saveFigure(figpath, self.figfile) | |
|
1168 | # store png plot to FTP server according to RT-Web format | |
|
1169 | name = self.getNameToFtp(thisDatetime, self.FTP_WEI, self.EXP_CODE, self.SUB_EXP_CODE, self.PLOT_CODE, self.PLOT_POS) | |
|
1170 | ftp_filename = os.path.join(figpath, name) | |
|
1171 | self.saveFigure(figpath, ftp_filename) | |
|
1172 | ||
|
1173 | self.counter_imagwr = 0 | |
|
1174 | ||
|
1175 | if x[1] >= self.axesList[0].xmax: | |
|
1176 | self.counter_imagwr = wr_period | |
|
1177 | self.__isConfig = False | |
|
774 | 1178 | self.figfile = None No newline at end of file |
@@ -316,7 +316,7 def createPmultilineYAxis(ax, x, y, xmin, xmax, ymin, ymax, xlabel='', ylabel='' | |||
|
316 | 316 | matplotlib.pyplot.ioff() |
|
317 | 317 | |
|
318 | 318 | # lines = ax.plot(x, y.T, marker=marker,markersize=markersize,linestyle=linestyle) |
|
319 |
lines = ax.plot(x, y.T, linestyle= |
|
|
319 | lines = ax.plot(x, y.T, linestyle=linestyle, marker=marker, markersize=markersize) | |
|
320 | 320 | leg = ax.legend(lines, legendlabels, loc='upper left', bbox_to_anchor=(1.01, 1.00), numpoints=1, handlelength=1.5, \ |
|
321 | 321 | handletextpad=0.5, borderpad=0.5, labelspacing=0.5, borderaxespad=0.) |
|
322 | 322 |
@@ -7,7 +7,9 from scipy import stats | |||
|
7 | 7 | import re |
|
8 | 8 | import datetime |
|
9 | 9 | import copy |
|
10 | ||
|
10 | import sys | |
|
11 | import importlib | |
|
12 | import itertools | |
|
11 | 13 | |
|
12 | 14 | from jroproc_base import ProcessingUnit, Operation |
|
13 | 15 | from model.data.jrodata import Parameters |
@@ -60,7 +62,7 class ParametersProc(ProcessingUnit): | |||
|
60 | 62 | |
|
61 | 63 | def run(self, nSeconds = None, nProfiles = None): |
|
62 | 64 | |
|
63 | self.dataOut.flagNoData = True | |
|
65 | ||
|
64 | 66 | |
|
65 | 67 | if self.firstdatatime == None: |
|
66 | 68 | self.firstdatatime = self.dataIn.utctime |
@@ -68,6 +70,7 class ParametersProc(ProcessingUnit): | |||
|
68 | 70 | #---------------------- Voltage Data --------------------------- |
|
69 | 71 | |
|
70 | 72 | if self.dataIn.type == "Voltage": |
|
73 | self.dataOut.flagNoData = True | |
|
71 | 74 | if nSeconds != None: |
|
72 | 75 | self.nSeconds = nSeconds |
|
73 | 76 | self.nProfiles= int(numpy.floor(nSeconds/(self.dataIn.ippSeconds*self.dataIn.nCohInt))) |
@@ -100,6 +103,7 class ParametersProc(ProcessingUnit): | |||
|
100 | 103 | self.dataOut.abscissaRange = self.dataIn.getVelRange(1) |
|
101 | 104 | self.dataOut.noise = self.dataIn.getNoise() |
|
102 | 105 | self.dataOut.normFactor = self.dataIn.normFactor |
|
106 | self.dataOut.flagNoData = False | |
|
103 | 107 | |
|
104 | 108 | #---------------------- Correlation Data --------------------------- |
|
105 | 109 | |
@@ -112,14 +116,14 class ParametersProc(ProcessingUnit): | |||
|
112 | 116 | self.dataOut.noise = self.dataIn.noise |
|
113 | 117 | self.dataOut.normFactor = self.dataIn.normFactor |
|
114 | 118 | self.dataOut.SNR = self.dataIn.SNR |
|
115 |
self.dataOut. |
|
|
119 | self.dataOut.groupList = self.dataIn.pairsList | |
|
120 | self.dataOut.flagNoData = False | |
|
116 | 121 | |
|
117 | 122 | |
|
118 | 123 | self.__updateObjFromInput() |
|
119 | self.dataOut.flagNoData = False | |
|
120 | 124 | self.firstdatatime = None |
|
121 | 125 | self.dataOut.initUtcTime = self.dataIn.ltctime |
|
122 |
self.dataOut. |
|
|
126 | self.dataOut.outputInterval = self.dataIn.timeInterval | |
|
123 | 127 | |
|
124 | 128 | #------------------- Get Moments ---------------------------------- |
|
125 | 129 | def GetMoments(self, channelList = None): |
@@ -238,7 +242,7 class ParametersProc(ProcessingUnit): | |||
|
238 | 242 | self.dataOut.noise |
|
239 | 243 | self.dataOut.normFactor |
|
240 | 244 | self.dataOut.SNR |
|
241 |
self.dataOut. |
|
|
245 | self.dataOut.groupList | |
|
242 | 246 | self.dataOut.nChannels |
|
243 | 247 | |
|
244 | 248 | Affected: |
@@ -251,7 +255,7 class ParametersProc(ProcessingUnit): | |||
|
251 | 255 | absc = self.dataOut.abscissaRange[:-1] |
|
252 | 256 | noise = self.dataOut.noise |
|
253 | 257 | SNR = self.dataOut.SNR |
|
254 |
pairsList = self.dataOut. |
|
|
258 | pairsList = self.dataOut.groupList | |
|
255 | 259 | nChannels = self.dataOut.nChannels |
|
256 | 260 | pairsAutoCorr, pairsCrossCorr = self.__getPairsAutoCorr(pairsList, nChannels) |
|
257 | 261 | self.dataOut.data_param = numpy.zeros((len(pairsCrossCorr)*2 + 1, nHeights)) |
@@ -1109,6 +1113,137 class ParametersProc(ProcessingUnit): | |||
|
1109 | 1113 | |
|
1110 | 1114 | return heights, error |
|
1111 | 1115 | |
|
1116 | def SpectralFitting(self, getSNR = True, path=None, file=None, groupList=None): | |
|
1117 | ||
|
1118 | ''' | |
|
1119 | Function GetMoments() | |
|
1120 | ||
|
1121 | Input: | |
|
1122 | Output: | |
|
1123 | Variables modified: | |
|
1124 | ''' | |
|
1125 | if path != None: | |
|
1126 | sys.path.append(path) | |
|
1127 | self.dataOut.library = importlib.import_module(file) | |
|
1128 | ||
|
1129 | #To be inserted as a parameter | |
|
1130 | groupArray = numpy.array(groupList) | |
|
1131 | # groupArray = numpy.array([[0,1],[2,3]]) | |
|
1132 | self.dataOut.groupList = groupArray | |
|
1133 | ||
|
1134 | nGroups = groupArray.shape[0] | |
|
1135 | nChannels = self.dataIn.nChannels | |
|
1136 | nHeights=self.dataIn.heightList.size | |
|
1137 | ||
|
1138 | #Parameters Array | |
|
1139 | self.dataOut.data_param = None | |
|
1140 | ||
|
1141 | #Set constants | |
|
1142 | constants = self.dataOut.library.setConstants(self.dataIn) | |
|
1143 | self.dataOut.constants = constants | |
|
1144 | M = self.dataIn.normFactor | |
|
1145 | N = self.dataIn.nFFTPoints | |
|
1146 | ippSeconds = self.dataIn.ippSeconds | |
|
1147 | K = self.dataIn.nIncohInt | |
|
1148 | pairsArray = numpy.array(self.dataIn.pairsList) | |
|
1149 | ||
|
1150 | #List of possible combinations | |
|
1151 | listComb = itertools.combinations(numpy.arange(groupArray.shape[1]),2) | |
|
1152 | indCross = numpy.zeros(len(list(listComb)), dtype = 'int') | |
|
1153 | ||
|
1154 | if getSNR: | |
|
1155 | listChannels = groupArray.reshape((groupArray.size)) | |
|
1156 | listChannels.sort() | |
|
1157 | noise = self.dataIn.getNoise() | |
|
1158 | self.dataOut.SNR = self.__getSNR(self.dataIn.data_spc[listChannels,:,:], noise[listChannels]) | |
|
1159 | ||
|
1160 | for i in range(nGroups): | |
|
1161 | coord = groupArray[i,:] | |
|
1162 | ||
|
1163 | #Input data array | |
|
1164 | data = self.dataIn.data_spc[coord,:,:]/(M*N) | |
|
1165 | data = data.reshape((data.shape[0]*data.shape[1],data.shape[2])) | |
|
1166 | ||
|
1167 | #Cross Spectra data array for Covariance Matrixes | |
|
1168 | ind = 0 | |
|
1169 | for pairs in listComb: | |
|
1170 | pairsSel = numpy.array([coord[x],coord[y]]) | |
|
1171 | indCross[ind] = int(numpy.where(numpy.all(pairsArray == pairsSel, axis = 1))[0][0]) | |
|
1172 | ind += 1 | |
|
1173 | dataCross = self.dataIn.data_cspc[indCross,:,:]/(M*N) | |
|
1174 | dataCross = dataCross**2/K | |
|
1175 | ||
|
1176 | for h in range(nHeights): | |
|
1177 | # print self.dataOut.heightList[h] | |
|
1178 | ||
|
1179 | #Input | |
|
1180 | d = data[:,h] | |
|
1181 | ||
|
1182 | #Covariance Matrix | |
|
1183 | D = numpy.diag(d**2/K) | |
|
1184 | ind = 0 | |
|
1185 | for pairs in listComb: | |
|
1186 | #Coordinates in Covariance Matrix | |
|
1187 | x = pairs[0] | |
|
1188 | y = pairs[1] | |
|
1189 | #Channel Index | |
|
1190 | S12 = dataCross[ind,:,h] | |
|
1191 | D12 = numpy.diag(S12) | |
|
1192 | #Completing Covariance Matrix with Cross Spectras | |
|
1193 | D[x*N:(x+1)*N,y*N:(y+1)*N] = D12 | |
|
1194 | D[y*N:(y+1)*N,x*N:(x+1)*N] = D12 | |
|
1195 | ind += 1 | |
|
1196 | Dinv=numpy.linalg.inv(D) | |
|
1197 | L=numpy.linalg.cholesky(Dinv) | |
|
1198 | LT=L.T | |
|
1199 | ||
|
1200 | dp = numpy.dot(LT,d) | |
|
1201 | ||
|
1202 | #Initial values | |
|
1203 | data_spc = self.dataIn.data_spc[coord,:,h] | |
|
1204 | p0 = self.dataOut.library.initialValuesFunction(data_spc, constants) | |
|
1205 | ||
|
1206 | #Least Squares | |
|
1207 | minp,covp,infodict,mesg,ier = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants),full_output=True) | |
|
1208 | # minp,covp = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants)) | |
|
1209 | #Chi square error | |
|
1210 | error0 = numpy.sum(infodict['fvec']**2)/(2*N) | |
|
1211 | # error0 = 0 | |
|
1212 | #Error with Jacobian | |
|
1213 | error1 = self.dataOut.library.errorFunction(minp,constants,LT) | |
|
1214 | #Save | |
|
1215 | if self.dataOut.data_param == None: | |
|
1216 | self.dataOut.data_param = numpy.zeros((nGroups, minp.size, nHeights))*numpy.nan | |
|
1217 | self.dataOut.error = numpy.zeros((nGroups, error1.size + 1, nHeights))*numpy.nan | |
|
1218 | ||
|
1219 | self.dataOut.error[i,:,h] = numpy.hstack((error0,error1)) | |
|
1220 | self.dataOut.data_param[i,:,h] = minp | |
|
1221 | return | |
|
1222 | ||
|
1223 | ||
|
1224 | def __residFunction(self, p, dp, LT, constants): | |
|
1225 | ||
|
1226 | fm = self.dataOut.library.modelFunction(p, constants) | |
|
1227 | fmp=numpy.dot(LT,fm) | |
|
1228 | ||
|
1229 | return dp-fmp | |
|
1230 | ||
|
1231 | def __getSNR(self, z, noise): | |
|
1232 | ||
|
1233 | avg = numpy.average(z, axis=1) | |
|
1234 | SNR = (avg.T-noise)/noise | |
|
1235 | SNR = SNR.T | |
|
1236 | return SNR | |
|
1237 | ||
|
1238 | def __chisq(p,chindex,hindex): | |
|
1239 | #similar to Resid but calculates CHI**2 | |
|
1240 | [LT,d,fm]=setupLTdfm(p,chindex,hindex) | |
|
1241 | dp=numpy.dot(LT,d) | |
|
1242 | fmp=numpy.dot(LT,fm) | |
|
1243 | chisq=numpy.dot((dp-fmp).T,(dp-fmp)) | |
|
1244 | return chisq | |
|
1245 | ||
|
1246 | ||
|
1112 | 1247 | |
|
1113 | 1248 | class WindProfiler(Operation): |
|
1114 | 1249 | |
@@ -1359,12 +1494,12 class WindProfiler(Operation): | |||
|
1359 | 1494 | winds = correctFactor*winds |
|
1360 | 1495 | return winds |
|
1361 | 1496 | |
|
1362 |
def __checkTime(self, currentTime, paramInterval, |
|
|
1497 | def __checkTime(self, currentTime, paramInterval, outputInterval): | |
|
1363 | 1498 | |
|
1364 | 1499 | dataTime = currentTime + paramInterval |
|
1365 | 1500 | deltaTime = dataTime - self.__initime |
|
1366 | 1501 | |
|
1367 |
if deltaTime >= |
|
|
1502 | if deltaTime >= outputInterval or deltaTime < 0: | |
|
1368 | 1503 | self.__dataReady = True |
|
1369 | 1504 | return |
|
1370 | 1505 | |
@@ -1462,7 +1597,7 class WindProfiler(Operation): | |||
|
1462 | 1597 | theta_y = theta_y[arrayChannel] |
|
1463 | 1598 | |
|
1464 | 1599 | velRadial0 = param[:,1,:] #Radial velocity |
|
1465 |
dataOut. |
|
|
1600 | dataOut.data_output, dataOut.heightRange, dataOut.SNR = self.techniqueDBS(velRadial0, theta_x, theta_y, azimuth, correctFactor, horizontalOnly, heightRange, SNR) #DBS Function | |
|
1466 | 1601 | |
|
1467 | 1602 | elif technique == 'SA': |
|
1468 | 1603 | |
@@ -1483,12 +1618,12 class WindProfiler(Operation): | |||
|
1483 | 1618 | |
|
1484 | 1619 | tau = dataOut.data_param |
|
1485 | 1620 | _lambda = dataOut.C/dataOut.frequency |
|
1486 |
pairsList = dataOut. |
|
|
1621 | pairsList = dataOut.groupList | |
|
1487 | 1622 | nChannels = dataOut.nChannels |
|
1488 | 1623 | |
|
1489 |
dataOut. |
|
|
1624 | dataOut.data_output = self.techniqueSA(pairs, pairsList, nChannels, tau, azimuth, _lambda, position_x, position_y, absc, correctFactor) | |
|
1490 | 1625 | dataOut.initUtcTime = dataOut.ltctime |
|
1491 |
dataOut. |
|
|
1626 | dataOut.outputInterval = dataOut.timeInterval | |
|
1492 | 1627 | |
|
1493 | 1628 | elif technique == 'Meteors': |
|
1494 | 1629 | dataOut.flagNoData = True |
@@ -1511,7 +1646,7 class WindProfiler(Operation): | |||
|
1511 | 1646 | hmax = kwargs['hmax'] |
|
1512 | 1647 | else: hmax = 110 |
|
1513 | 1648 | |
|
1514 |
dataOut. |
|
|
1649 | dataOut.outputInterval = nHours*3600 | |
|
1515 | 1650 | |
|
1516 | 1651 | if self.__isConfig == False: |
|
1517 | 1652 | # self.__initime = dataOut.datatime.replace(minute = 0, second = 0, microsecond = 03) |
@@ -1526,14 +1661,89 class WindProfiler(Operation): | |||
|
1526 | 1661 | else: |
|
1527 | 1662 | self.__buffer = numpy.vstack((self.__buffer, dataOut.data_param)) |
|
1528 | 1663 | |
|
1529 |
self.__checkTime(dataOut.ltctime, dataOut.paramInterval, dataOut. |
|
|
1664 | self.__checkTime(dataOut.ltctime, dataOut.paramInterval, dataOut.outputInterval) #Check if the buffer is ready | |
|
1530 | 1665 | |
|
1531 | 1666 | if self.__dataReady: |
|
1532 | 1667 | dataOut.initUtcTime = self.__initime |
|
1533 |
self.__initime = self.__initime + dataOut. |
|
|
1668 | self.__initime = self.__initime + dataOut.outputInterval #to erase time offset | |
|
1534 | 1669 | |
|
1535 |
dataOut. |
|
|
1670 | dataOut.data_output, dataOut.heightRange = self.techniqueMeteors(self.__buffer, meteorThresh, hmin, hmax) | |
|
1536 | 1671 | dataOut.flagNoData = False |
|
1537 | 1672 | self.__buffer = None |
|
1538 | 1673 | |
|
1539 | return No newline at end of file | |
|
1674 | return | |
|
1675 | ||
|
1676 | class EWDriftsEstimation(Operation): | |
|
1677 | ||
|
1678 | ||
|
1679 | def __init__(self): | |
|
1680 | Operation.__init__(self) | |
|
1681 | ||
|
1682 | def __correctValues(self, heiRang, phi, velRadial, SNR): | |
|
1683 | listPhi = phi.tolist() | |
|
1684 | maxid = listPhi.index(max(listPhi)) | |
|
1685 | minid = listPhi.index(min(listPhi)) | |
|
1686 | ||
|
1687 | rango = range(len(phi)) | |
|
1688 | # rango = numpy.delete(rango,maxid) | |
|
1689 | ||
|
1690 | heiRang1 = heiRang*math.cos(phi[maxid]) | |
|
1691 | heiRangAux = heiRang*math.cos(phi[minid]) | |
|
1692 | indOut = (heiRang1 < heiRangAux[0]).nonzero() | |
|
1693 | heiRang1 = numpy.delete(heiRang1,indOut) | |
|
1694 | ||
|
1695 | velRadial1 = numpy.zeros([len(phi),len(heiRang1)]) | |
|
1696 | SNR1 = numpy.zeros([len(phi),len(heiRang1)]) | |
|
1697 | ||
|
1698 | for i in rango: | |
|
1699 | x = heiRang*math.cos(phi[i]) | |
|
1700 | y1 = velRadial[i,:] | |
|
1701 | f1 = interpolate.interp1d(x,y1,kind = 'cubic') | |
|
1702 | ||
|
1703 | x1 = heiRang1 | |
|
1704 | y11 = f1(x1) | |
|
1705 | ||
|
1706 | y2 = SNR[i,:] | |
|
1707 | f2 = interpolate.interp1d(x,y2,kind = 'cubic') | |
|
1708 | y21 = f2(x1) | |
|
1709 | ||
|
1710 | velRadial1[i,:] = y11 | |
|
1711 | SNR1[i,:] = y21 | |
|
1712 | ||
|
1713 | return heiRang1, velRadial1, SNR1 | |
|
1714 | ||
|
1715 | def run(self, dataOut, zenith, zenithCorrection): | |
|
1716 | heiRang = dataOut.heightList | |
|
1717 | velRadial = dataOut.data_param[:,3,:] | |
|
1718 | SNR = dataOut.SNR | |
|
1719 | ||
|
1720 | zenith = numpy.array(zenith) | |
|
1721 | zenith -= zenithCorrection | |
|
1722 | zenith *= numpy.pi/180 | |
|
1723 | ||
|
1724 | heiRang1, velRadial1, SNR1 = self.__correctValues(heiRang, numpy.abs(zenith), velRadial, SNR) | |
|
1725 | ||
|
1726 | alp = zenith[0] | |
|
1727 | bet = zenith[1] | |
|
1728 | ||
|
1729 | w_w = velRadial1[0,:] | |
|
1730 | w_e = velRadial1[1,:] | |
|
1731 | ||
|
1732 | w = (w_w*numpy.sin(bet) - w_e*numpy.sin(alp))/(numpy.cos(alp)*numpy.sin(bet) - numpy.cos(bet)*numpy.sin(alp)) | |
|
1733 | u = (w_w*numpy.cos(bet) - w_e*numpy.cos(alp))/(numpy.sin(alp)*numpy.cos(bet) - numpy.sin(bet)*numpy.cos(alp)) | |
|
1734 | ||
|
1735 | winds = numpy.vstack((u,w)) | |
|
1736 | ||
|
1737 | dataOut.heightList = heiRang1 | |
|
1738 | dataOut.data_output = winds | |
|
1739 | dataOut.SNR = SNR1 | |
|
1740 | ||
|
1741 | dataOut.initUtcTime = dataOut.ltctime | |
|
1742 | dataOut.outputInterval = dataOut.timeInterval | |
|
1743 | return | |
|
1744 | ||
|
1745 | ||
|
1746 | ||
|
1747 | ||
|
1748 | ||
|
1749 | No newline at end of file |
General Comments 0
You need to be logged in to leave comments.
Login now