##// END OF EJS Templates
schain
RSS
Drifts proccesing for some heights of interest. Some parameters are needded for madrigal writting: lat, lon. NAN vectors are not written to avoid errors at madrigal view.
imanay -
r1690:571f81d70cd0
parent child
Show More
Side by Side Unified
Context file:
r1690:571f81d70cd0 Loading diff...:
Show file before | Show file after | Hide comments
@@ -1455,6 +1455,7 class SpectralMoments(Operation):
1455
1455
1456 for ind in range(nChannel):
1456 for ind in range(nChannel):
1457 data_param[ind,:,:] = self.__calculateMoments( data[ind,:,:] , absc , noise[ind], nicoh=nIncohInt, smooth=smooth, type1=type1, fwindow=fwindow)
1457 data_param[ind,:,:] = self.__calculateMoments( data[ind,:,:] , absc , noise[ind], nicoh=nIncohInt, smooth=smooth, type1=type1, fwindow=fwindow)
1458 #print('snr:',data_param[:,0])
1458
1459
1459 if proc_type == 1:
1460 if proc_type == 1:
1460 dataOut.moments = data_param[:,1:,:]
1461 dataOut.moments = data_param[:,1:,:]
@@ -1963,68 +1964,52 class JULIADriftsEstimation(Operation):
1963
1964
1964 def run(self, dataOut, zenith, zenithCorrection,heights=None, statistics=0, otype=0):
1965 def run(self, dataOut, zenith, zenithCorrection,heights=None, statistics=0, otype=0):
1965
1966
1966 nCh=dataOut.spcpar.shape[0]
1967
1968 dataOut.lat=-11.95
1969 dataOut.lon=-76.87
1967
1970
1971 nCh=dataOut.spcpar.shape[0]
1968 nHei=dataOut.spcpar.shape[1]
1972 nHei=dataOut.spcpar.shape[1]
1969 nParam=dataOut.spcpar.shape[2]
1973 nParam=dataOut.spcpar.shape[2]
1970 # Solo las alturas de interes
1974 # SelecciΓ³n de alturas
1971 hei=dataOut.heightList
1975
1972 hvalid=numpy.where([hei >= heights[0]][0] & [hei <= heights[1]][0])[0]
1976 if not heights:
1973 nhvalid=len(hvalid)
1977 parm = numpy.zeros((nCh,nHei,nParam))
1974 parm = numpy.zeros((nCh,nhvalid,nParam))
1978 parm[:] = dataOut.spcpar[:]
1975 parm = dataOut.spcpar[:,hvalid,:]
1979 else:
1980 hei=dataOut.heightList
1981 hvalid=numpy.where([hei >= heights[0]][0] & [hei <= heights[1]][0])[0]
1982 nhvalid=len(hvalid)
1983 dataOut.heightList = hei[hvalid]
1984 parm = numpy.zeros((nCh,nhvalid,nParam))
1985 parm[:] = dataOut.spcpar[:,hvalid,:]
1986
1976
1987
1977 # Primer filtrado: Umbral de SNR
1988 # Primer filtrado: Umbral de SNR
1978 #snrth=-19
1979 for i in range(nCh):
1989 for i in range(nCh):
1980 #print('snr:',parm[i,:,2])
1990 parm[i,:,:] = self.data_filter(parm[i,:,:])[0]
1981 #dataOut.spcpar[i,hvalid,:] = self.data_filter(parm[i,:,:],snrth)[0]
1991
1982 dataOut.spcpar[i,hvalid,:] = self.data_filter(parm[i,:,:])[0]
1983 #print('dataOut.spcpar[0,:,2]',dataOut.spcpar[0,:,2])
1984 #print('dataOut.spcpar[1,:,2]',dataOut.spcpar[1,:,2])
1985 zenith = numpy.array(zenith)
1992 zenith = numpy.array(zenith)
1986 zenith -= zenithCorrection
1993 zenith -= zenithCorrection
1987 zenith *= numpy.pi/180
1994 zenith *= numpy.pi/180
1988 alpha = zenith[0]
1995 alpha = zenith[0]
1989 beta = zenith[1]
1996 beta = zenith[1]
1990
1997 dopplerCH0 = parm[0,:,0]
1991 dopplerCH0 = dataOut.spcpar[0,:,0]
1998 dopplerCH1 = parm[1,:,0]
1992 dopplerCH1 = dataOut.spcpar[1,:,0]
1999 swCH0 = parm[0,:,1]
1993 swCH0 = dataOut.spcpar[0,:,1]
2000 swCH1 = parm[1,:,1]
1994 swCH1 = dataOut.spcpar[1,:,1]
2001 snrCH0 = 10*numpy.log10(parm[0,:,2])
1995 snrCH0 = 10*numpy.log10(dataOut.spcpar[0,:,2])
2002 snrCH1 = 10*numpy.log10(parm[1,:,2])
1996 snrCH1 = 10*numpy.log10(dataOut.spcpar[1,:,2])
2003 noiseCH0 = parm[0,:,3]
1997 noiseCH0 = dataOut.spcpar[0,:,3]
2004 noiseCH1 = parm[1,:,3]
1998 noiseCH1 = dataOut.spcpar[1,:,3]
2005 wErrCH0 = parm[0,:,5]
1999 wErrCH0 = dataOut.spcpar[0,:,5]
2006 wErrCH1 = parm[1,:,5]
2000 wErrCH1 = dataOut.spcpar[1,:,5]
2001
2007
2002 # Vertical and zonal calculation according to geometry
2008 # Vertical and zonal calculation according to geometry
2003 sinB_A = numpy.sin(beta)*numpy.cos(alpha) - numpy.sin(alpha)* numpy.cos(beta)
2009 sinB_A = numpy.sin(beta)*numpy.cos(alpha) - numpy.sin(alpha)* numpy.cos(beta)
2004 drift = -(dopplerCH0 * numpy.sin(beta) - dopplerCH1 * numpy.sin(alpha))/ sinB_A
2010 drift = -(dopplerCH0 * numpy.sin(beta) - dopplerCH1 * numpy.sin(alpha))/ sinB_A
2005 '''
2006 print('drift.shape:',drift.shape)
2007 print('drift min:', numpy.nanmin(drift))
2008 print('drift max:', numpy.nanmax(drift))
2009 '''
2010
2011 '''
2012 print('shape:', dopplerCH0[hvalid].shape)
2013 print('dopplerCH0:', dopplerCH0[hvalid])
2014 print('dopplerCH1:', dopplerCH1[hvalid])
2015 print('drift:', drift[hvalid])
2016 '''
2017 zonal = (dopplerCH0 * numpy.cos(beta) - dopplerCH1 * numpy.cos(alpha))/ sinB_A
2011 zonal = (dopplerCH0 * numpy.cos(beta) - dopplerCH1 * numpy.cos(alpha))/ sinB_A
2018 '''
2019 print('zonal min:', numpy.nanmin(zonal))
2020 print('zonal max:', numpy.nanmax(zonal))
2021 '''
2022 #print('zonal:', zonal[hvalid])
2023 snr = (snrCH0 + snrCH1)/2
2012 snr = (snrCH0 + snrCH1)/2
2024 '''
2025 print('snr min:', 10*numpy.log10(numpy.nanmin(snr)))
2026 print('snr max:', 10*numpy.log10(numpy.nanmax(snr)))
2027 '''
2028 noise = (noiseCH0 + noiseCH1)/2
2013 noise = (noiseCH0 + noiseCH1)/2
2029 sw = (swCH0 + swCH1)/2
2014 sw = (swCH0 + swCH1)/2
2030 w_w_err= numpy.sqrt(numpy.power(wErrCH0 * numpy.sin(beta)/numpy.abs(sinB_A),2) + numpy.power(wErrCH1 * numpy.sin(alpha)/numpy.abs(sinB_A),2))
2015 w_w_err= numpy.sqrt(numpy.power(wErrCH0 * numpy.sin(beta)/numpy.abs(sinB_A),2) + numpy.power(wErrCH1 * numpy.sin(alpha)/numpy.abs(sinB_A),2))
@@ -2033,7 +2018,7 class JULIADriftsEstimation(Operation):
2033 # for statistics150km
2018 # for statistics150km
2034 if statistics:
2019 if statistics:
2035 print('Implemented offline.')
2020 print('Implemented offline.')
2036
2021
2037 if otype == 0:
2022 if otype == 0:
2038 winds = numpy.vstack((snr, drift, zonal, noise, sw, w_w_err, w_e_err)) # to process statistics drifts
2023 winds = numpy.vstack((snr, drift, zonal, noise, sw, w_w_err, w_e_err)) # to process statistics drifts
2039 elif otype == 3:
2024 elif otype == 3:
@@ -2042,13 +2027,14 class JULIADriftsEstimation(Operation):
2042 winds = numpy.vstack((snrCH0, drift, snrCH1, zonal)) # to generic plot: 4 RTI's
2027 winds = numpy.vstack((snrCH0, drift, snrCH1, zonal)) # to generic plot: 4 RTI's
2043
2028
2044 snr1 = numpy.vstack((snrCH0, snrCH1))
2029 snr1 = numpy.vstack((snrCH0, snrCH1))
2045
2046 dataOut.data_output = winds
2030 dataOut.data_output = winds
2047 dataOut.data_snr = snr1
2031 dataOut.data_snr = snr1
2048
2032
2049 dataOut.utctimeInit = dataOut.utctime
2033 dataOut.utctimeInit = dataOut.utctime
2050 dataOut.outputInterval = dataOut.timeInterval
2034 dataOut.outputInterval = dataOut.timeInterval
2051
2035
2036 dataOut.flagNoData = numpy.all(numpy.isnan(dataOut.data_output[0])) # NAN vectors are not written
2037
2052 return dataOut
2038 return dataOut
2053
2039
2054 class SALags(Operation):
2040 class SALags(Operation):
General Comments 0
You need to be logged in to leave comments. Login now