##// END OF EJS Templates
schain
RSS
test
Alexander Valdez -
r1681:82935db89bc7
parent child
Show More
Side by Side Unified
Context file:
r1681:82935db89bc7 Loading diff...:
Show file before | Show file after | Hide comments
@@ -2855,16 +2855,22 class SpectralFitting(Operation):
2855 p0 = dataOut.data_param[i,:,h-1]
2855 p0 = dataOut.data_param[i,:,h-1]
2856 else:
2856 else:
2857 p0 = numpy.array(self.library.initialValuesFunction(data_spc*w, constants))# sin el i(data_spc, constants, i)
2857 p0 = numpy.array(self.library.initialValuesFunction(data_spc*w, constants))# sin el i(data_spc, constants, i)
2858 #print("AQUI REEMPLAZAS CON EL fd0",fd0)
2858 p0[3] = fd0
2859 p0[3] = fd0
2859 if filec != None:
2860 if filec != None:
2860 p0 = self.weightf.Vrfit(p0,tini.tm_year,tini.tm_yday,index,h,i)
2861 p0 = self.weightf.Vrfit(p0,tini.tm_year,tini.tm_yday,index,h,i)
2862
2863 #print("minP0-ANTES DE OPTIMIZE")
2864 #print(p0)
2861 try:
2865 try:
2862 #Least Squares
2866 #Least Squares
2863 minp,covp,infodict,mesg,ier = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants),full_output=True)
2867 minp,covp,infodict,mesg,ier = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants),full_output=True)
2864 #minp,covp = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants))
2868 print("VERIFICAR SI ESTA VARIANDO minp--------------------------------", minp[3])
2865 #Chi square error
2869 print("COMPARACION ---- p0[3],fd0",p0[3],fd0)
2870 #minp,covp = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants))
2871 #Chi square error
2866 error0 = numpy.sum(infodict['fvec']**2)/(2*N)
2872 error0 = numpy.sum(infodict['fvec']**2)/(2*N)
2867 #Error with Jacobian
2873 #Error with Jacobian
2868 error1 = self.library.errorFunction(minp,constants,LT)
2874 error1 = self.library.errorFunction(minp,constants,LT)
2869
2875
2870 except:
2876 except:
@@ -2880,7 +2886,7 class SpectralFitting(Operation):
2880 if dataOut.data_param is None:
2886 if dataOut.data_param is None:
2881 dataOut.data_param = numpy.zeros((nGroups, p0.size, nHeights))*numpy.nan
2887 dataOut.data_param = numpy.zeros((nGroups, p0.size, nHeights))*numpy.nan
2882 dataOut.data_error = numpy.zeros((nGroups, p0.size + 1, nHeights))*numpy.nan
2888 dataOut.data_error = numpy.zeros((nGroups, p0.size + 1, nHeights))*numpy.nan
2883
2889 print("CLASE SF minp- ?????",minp)
2884 dataOut.data_error[i,:,h] = numpy.hstack((error0,error1))
2890 dataOut.data_error[i,:,h] = numpy.hstack((error0,error1))
2885 dataOut.data_param[i,:,h] = minp
2891 dataOut.data_param[i,:,h] = minp
2886 for ht in range(nHeights-1) :
2892 for ht in range(nHeights-1) :
@@ -3648,6 +3654,7 class EWDriftsEstimation(Operation):
3648 else :
3654 else :
3649 moments=numpy.vstack(([velRadialm[0,:]],[velRadialm[0,:]]))
3655 moments=numpy.vstack(([velRadialm[0,:]],[velRadialm[0,:]]))
3650 dataOut.moments=moments
3656 dataOut.moments=moments
3657 print("CLASE EWD Estimation",moments)
3651 # Coherent
3658 # Coherent
3652 smooth_wC = ebufc[0,:]
3659 smooth_wC = ebufc[0,:]
3653 p_w0C = rbufc[0,:]
3660 p_w0C = rbufc[0,:]
@@ -5034,127 +5041,10 class SMOperations():
5034 chan1Y = chan1
5041 chan1Y = chan1
5035 chan2Y = chan2
5042 chan2Y = chan2
5036 distances[2:4] = numpy.array([d1,d2])
5043 distances[2:4] = numpy.array([d1,d2])
5037 # axisXsides = numpy.reshape(axisX[ix,:],4)
5044
5038 #
5039 # channelCentX = int(numpy.intersect1d(pairX[0,:], pairX[1,:])[0])
5040 # channelCentY = int(numpy.intersect1d(pairY[0,:], pairY[1,:])[0])
5041 #
5042 # ind25X = numpy.where(pairX[0,:] != channelCentX)[0][0]
5043 # ind20X = numpy.where(pairX[1,:] != channelCentX)[0][0]
5044 # channel25X = int(pairX[0,ind25X])
5045 # channel20X = int(pairX[1,ind20X])
5046 # ind25Y = numpy.where(pairY[0,:] != channelCentY)[0][0]
5047 # ind20Y = numpy.where(pairY[1,:] != channelCentY)[0][0]
5048 # channel25Y = int(pairY[0,ind25Y])
5049 # channel20Y = int(pairY[1,ind20Y])
5050
5051 # pairslist = [(channelCentX, channel25X),(channelCentX, channel20X),(channelCentY,channel25Y),(channelCentY, channel20Y)]
5052 pairslist = [(chanCX, chan1X),(chanCX, chan2X),(chanCY,chan1Y),(chanCY, chan2Y)]
5045 pairslist = [(chanCX, chan1X),(chanCX, chan2X),(chanCY,chan1Y),(chanCY, chan2Y)]
5053
5046
5054 return pairslist, distances
5047 return pairslist, distances
5055 # def __getAOA(self, phases, pairsList, error, AOAthresh, azimuth):
5056 #
5057 # arrayAOA = numpy.zeros((phases.shape[0],3))
5058 # cosdir0, cosdir = self.__getDirectionCosines(phases, pairsList)
5059 #
5060 # arrayAOA[:,:2] = self.__calculateAOA(cosdir, azimuth)
5061 # cosDirError = numpy.sum(numpy.abs(cosdir0 - cosdir), axis = 1)
5062 # arrayAOA[:,2] = cosDirError
5063 #
5064 # azimuthAngle = arrayAOA[:,0]
5065 # zenithAngle = arrayAOA[:,1]
5066 #
5067 # #Setting Error
5068 # #Number 3: AOA not fesible
5069 # indInvalid = numpy.where(numpy.logical_and((numpy.logical_or(numpy.isnan(zenithAngle), numpy.isnan(azimuthAngle))),error == 0))[0]
5070 # error[indInvalid] = 3
5071 # #Number 4: Large difference in AOAs obtained from different antenna baselines
5072 # indInvalid = numpy.where(numpy.logical_and(cosDirError > AOAthresh,error == 0))[0]
5073 # error[indInvalid] = 4
5074 # return arrayAOA, error
5075 #
5076 # def __getDirectionCosines(self, arrayPhase, pairsList):
5077 #
5078 # #Initializing some variables
5079 # ang_aux = numpy.array([-8,-7,-6,-5,-4,-3,-2,-1,0,1,2,3,4,5,6,7,8])*2*numpy.pi
5080 # ang_aux = ang_aux.reshape(1,ang_aux.size)
5081 #
5082 # cosdir = numpy.zeros((arrayPhase.shape[0],2))
5083 # cosdir0 = numpy.zeros((arrayPhase.shape[0],2))
5084 #
5085 #
5086 # for i in range(2):
5087 # #First Estimation
5088 # phi0_aux = arrayPhase[:,pairsList[i][0]] + arrayPhase[:,pairsList[i][1]]
5089 # #Dealias
5090 # indcsi = numpy.where(phi0_aux > numpy.pi)
5091 # phi0_aux[indcsi] -= 2*numpy.pi
5092 # indcsi = numpy.where(phi0_aux < -numpy.pi)
5093 # phi0_aux[indcsi] += 2*numpy.pi
5094 # #Direction Cosine 0
5095 # cosdir0[:,i] = -(phi0_aux)/(2*numpy.pi*0.5)
5096 #
5097 # #Most-Accurate Second Estimation
5098 # phi1_aux = arrayPhase[:,pairsList[i][0]] - arrayPhase[:,pairsList[i][1]]
5099 # phi1_aux = phi1_aux.reshape(phi1_aux.size,1)
5100 # #Direction Cosine 1
5101 # cosdir1 = -(phi1_aux + ang_aux)/(2*numpy.pi*4.5)
5102 #
5103 # #Searching the correct Direction Cosine
5104 # cosdir0_aux = cosdir0[:,i]
5105 # cosdir0_aux = cosdir0_aux.reshape(cosdir0_aux.size,1)
5106 # #Minimum Distance
5107 # cosDiff = (cosdir1 - cosdir0_aux)**2
5108 # indcos = cosDiff.argmin(axis = 1)
5109 # #Saving Value obtained
5110 # cosdir[:,i] = cosdir1[numpy.arange(len(indcos)),indcos]
5111 #
5112 # return cosdir0, cosdir
5113 #
5114 # def __calculateAOA(self, cosdir, azimuth):
5115 # cosdirX = cosdir[:,0]
5116 # cosdirY = cosdir[:,1]
5117 #
5118 # zenithAngle = numpy.arccos(numpy.sqrt(1 - cosdirX**2 - cosdirY**2))*180/numpy.pi
5119 # azimuthAngle = numpy.arctan2(cosdirX,cosdirY)*180/numpy.pi + azimuth #0 deg north, 90 deg east
5120 # angles = numpy.vstack((azimuthAngle, zenithAngle)).transpose()
5121 #
5122 # return angles
5123 #
5124 # def __getHeights(self, Ranges, zenith, error, minHeight, maxHeight):
5125 #
5126 # Ramb = 375 #Ramb = c/(2*PRF)
5127 # Re = 6371 #Earth Radius
5128 # heights = numpy.zeros(Ranges.shape)
5129 #
5130 # R_aux = numpy.array([0,1,2])*Ramb
5131 # R_aux = R_aux.reshape(1,R_aux.size)
5132 #
5133 # Ranges = Ranges.reshape(Ranges.size,1)
5134 #
5135 # Ri = Ranges + R_aux
5136 # hi = numpy.sqrt(Re**2 + Ri**2 + (2*Re*numpy.cos(zenith*numpy.pi/180)*Ri.transpose()).transpose()) - Re
5137 #
5138 # #Check if there is a height between 70 and 110 km
5139 # h_bool = numpy.sum(numpy.logical_and(hi > minHeight, hi < maxHeight), axis = 1)
5140 # ind_h = numpy.where(h_bool == 1)[0]
5141 #
5142 # hCorr = hi[ind_h, :]
5143 # ind_hCorr = numpy.where(numpy.logical_and(hi > minHeight, hi < maxHeight))
5144 #
5145 # hCorr = hi[ind_hCorr]
5146 # heights[ind_h] = hCorr
5147 #
5148 # #Setting Error
5149 # #Number 13: Height unresolvable echo: not valid height within 70 to 110 km
5150 # #Number 14: Height ambiguous echo: more than one possible height within 70 to 110 km
5151 #
5152 # indInvalid2 = numpy.where(numpy.logical_and(h_bool > 1, error == 0))[0]
5153 # error[indInvalid2] = 14
5154 # indInvalid1 = numpy.where(numpy.logical_and(h_bool == 0, error == 0))[0]
5155 # error[indInvalid1] = 13
5156 #
5157 # return heights, error
5158
5048
5159 class IGRFModel(Operation):
5049 class IGRFModel(Operation):
5160 """Operation to calculate Geomagnetic parameters.
5050 """Operation to calculate Geomagnetic parameters.
General Comments 0
You need to be logged in to leave comments. Login now