##// END OF EJS Templates
schain
RSS
LAST UPDATE AFTER PERCY & ALEX REVISION
Alexander Valdez -
r1672:65c1739783d6
parent child
Show More
Side by Side Unified
Context file:
r1672:65c1739783d6 Loading diff...:
Show file before | Show file after | Hide comments
This diff has been collapsed as it changes many lines, (602 lines changed) Show them Hide them
@@ -85,10 +85,11 class ParametersProc(ProcessingUnit):
85 85 self.dataOut.timeInterval1 = self.dataIn.timeInterval
86 86 self.dataOut.heightList = self.dataIn.heightList
87 87 self.dataOut.frequency = self.dataIn.frequency
88 self.dataOut.runNextUnit = self.dataIn.runNextUnit
88 89
90 def run(self, runNextUnit=0):
89 91
90 def run(self):
91
92 self.dataIn.runNextUnit = runNextUnit
92 93 #---------------------- Voltage Data ---------------------------
93 94
94 95 if self.dataIn.type == "Voltage":
@@ -882,7 +883,7 class PrecipitationProc(Operation):
882 883 self.i=0
883 884
884 885 def run(self, dataOut, radar=None, Pt=5000, Gt=295.1209, Gr=70.7945, Lambda=0.6741, aL=2.5118,
885 tauW=4e-06, ThetaT=0.1656317, ThetaR=0.36774087, Km2 = 0.93, Altitude=3350, SNRdBlimit=-30,
886 tauW=4e-06, ThetaT=0.1656317, ThetaR=0.36774087, Km2 = 0.93, Altitude=3350,SNRdBlimit=-30,
886 887 channel=None):
887 888
888 889 # print ('Entering PrecepitationProc ... ')
@@ -1454,7 +1455,6 class SpectralMoments(Operation):
1454 1455 dataOut.data_snr = data_param[:,0]
1455 1456 dataOut.data_pow = data_param[:,6] # to compare with type0 proccessing
1456 1457 dataOut.spcpar=numpy.stack((dataOut.data_dop,dataOut.data_width,dataOut.data_snr, data_param[:,3], data_param[:,4],data_param[:,5]),axis=2)
1457
1458 1458 else:
1459 1459 dataOut.moments = data_param[:,1:,:]
1460 1460 dataOut.data_snr = data_param[:,0]
@@ -2223,6 +2223,7 class SpectralFitting(Operation):
2223 2223 return moments
2224 2224
2225 2225 def __DiffCoherent(self, spectra, cspectra, dataOut, noise, snrth, coh_th, hei_th):
2226
2226 2227 nProf = dataOut.nProfiles
2227 2228 heights = dataOut.heightList
2228 2229 nHei = len(heights)
@@ -2248,20 +2249,17 class SpectralFitting(Operation):
2248 2249
2249 2250 if coh_th == None : coh_th = numpy.array([0.75,0.65,0.15]) # 0.65
2250 2251 if hei_th == None : hei_th = numpy.array([60,300,650])
2251 for ic in range(2):
2252 for ic in range(nPairs):
2252 2253 pair = crosspairs[ic]
2253 2254 #si el SNR es mayor que el SNR threshold los datos se toman coherentes
2254 2255 s_n0 = power[pair[0],:]/noise[pair[0]]
2255 2256 s_n1 = power[pair[1],:]/noise[pair[1]]
2256
2257 2257 valid1 =(s_n0>=snr_th).nonzero()
2258 2258 valid2 = (s_n1>=snr_th).nonzero()
2259 #valid = valid2 + valid1 #numpy.concatenate((valid1,valid2), axis=None)
2260 2259 valid1 = numpy.array(valid1[0])
2261 2260 valid2 = numpy.array(valid2[0])
2262 2261 valid = valid1
2263 2262 for iv in range(len(valid2)):
2264 #for ivv in range(len(valid1)) :
2265 2263 indv = numpy.array((valid1 == valid2[iv]).nonzero())
2266 2264 if len(indv[0]) == 0 :
2267 2265 valid = numpy.concatenate((valid,valid2[iv]), axis=None)
@@ -2269,17 +2267,13 class SpectralFitting(Operation):
2269 2267 my_coh_aver[pair[0],valid]=1
2270 2268 my_coh_aver[pair[1],valid]=1
2271 2269 # si la coherencia es mayor a la coherencia threshold los datos se toman
2272 #print my_coh_aver[0,:]
2273 2270 coh = numpy.squeeze(numpy.nansum(cspectra[ic,:,:], axis=0)/numpy.sqrt(numpy.nansum(spectra[pair[0],:,:], axis=0)*numpy.nansum(spectra[pair[1],:,:], axis=0)))
2274 #print('coh',numpy.absolute(coh))
2275 2271 for ih in range(len(hei_th)):
2276 2272 hvalid = (heights>hei_th[ih]).nonzero()
2277 2273 hvalid = hvalid[0]
2278 2274 if len(hvalid)>0:
2279 2275 valid = (numpy.absolute(coh[hvalid])>coh_th[ih]).nonzero()
2280 2276 valid = valid[0]
2281 #print('hvalid:',hvalid)
2282 #print('valid', valid)
2283 2277 if len(valid)>0:
2284 2278 my_coh_aver[pair[0],hvalid[valid]] =1
2285 2279 my_coh_aver[pair[1],hvalid[valid]] =1
@@ -2295,7 +2289,7 class SpectralFitting(Operation):
2295 2289 my_incoh_aver[pair[1],incoh_echoes] = 1
2296 2290
2297 2291
2298 for ic in range(2):
2292 for ic in range(nPairs):
2299 2293 pair = crosspairs[ic]
2300 2294
2301 2295 valid1 =(my_coh_aver[pair[0],:]==1 ).nonzero()
@@ -2303,29 +2297,24 class SpectralFitting(Operation):
2303 2297 valid1 = numpy.array(valid1[0])
2304 2298 valid2 = numpy.array(valid2[0])
2305 2299 valid = valid1
2306 #print valid1 , valid2
2300
2307 2301 for iv in range(len(valid2)):
2308 #for ivv in range(len(valid1)) :
2302
2309 2303 indv = numpy.array((valid1 == valid2[iv]).nonzero())
2310 2304 if len(indv[0]) == 0 :
2311 2305 valid = numpy.concatenate((valid,valid2[iv]), axis=None)
2312 #print valid
2313 #valid = numpy.concatenate((valid1,valid2), axis=None)
2314 2306 valid1 =(my_coh_aver[pair[0],:] !=1 ).nonzero()
2315 2307 valid2 = (my_coh_aver[pair[1],:] !=1).nonzero()
2316 2308 valid1 = numpy.array(valid1[0])
2317 2309 valid2 = numpy.array(valid2[0])
2318 2310 incoh_echoes = valid1
2319 #print valid1, valid2
2320 #incoh_echoes= numpy.concatenate((valid1,valid2), axis=None)
2321 2311 for iv in range(len(valid2)):
2322 #for ivv in range(len(valid1)) :
2312
2323 2313 indv = numpy.array((valid1 == valid2[iv]).nonzero())
2324 2314 if len(indv[0]) == 0 :
2325 2315 incoh_echoes = numpy.concatenate(( incoh_echoes,valid2[iv]), axis=None)
2326 #print incoh_echoes
2316
2327 2317 if len(valid)>0:
2328 #print pair
2329 2318 coh_spectra[pair[0],:,valid] = spectra[pair[0],:,valid]
2330 2319 coh_spectra[pair[1],:,valid] = spectra[pair[1],:,valid]
2331 2320 coh_cspectra[ic,:,valid] = cspectra[ic,:,valid]
@@ -2339,6 +2328,7 class SpectralFitting(Operation):
2339 2328 incoh_aver[pair[1],incoh_echoes]=1
2340 2329 return my_incoh_spectra ,my_incoh_cspectra,my_incoh_aver,my_coh_aver, incoh_spectra, coh_spectra, incoh_cspectra, coh_cspectra, incoh_aver, coh_aver
2341 2330
2331
2342 2332 def __CleanCoherent(self,snrth, spectra, cspectra, coh_aver,dataOut, noise,clean_coh_echoes,index):
2343 2333
2344 2334 nProf = dataOut.nProfiles
@@ -2349,10 +2339,7 class SpectralFitting(Operation):
2349 2339 crosspairs = dataOut.groupList
2350 2340 nPairs = len(crosspairs)
2351 2341
2352 #data = dataOut.data_pre[0]
2353 2342 absc = dataOut.abscissaList[:-1]
2354 #noise = dataOut.noise
2355 #nChannel = data.shape[0]
2356 2343 data_param = numpy.zeros((nChan, 4, spectra.shape[2]))
2357 2344 clean_coh_spectra = spectra.copy()
2358 2345 clean_coh_cspectra = cspectra.copy()
@@ -2364,17 +2351,12 class SpectralFitting(Operation):
2364 2351 rtime0 = [6,18] # periodo sin ESF
2365 2352 rtime1 = [10.5,13.5] # periodo con alta coherencia y alto ancho espectral (esperado): SOL.
2366 2353
2367 time = index*5./60
2354 time = index*5./60 # en base a 5 min de proceso
2368 2355 if clean_coh_echoes == 1 :
2369 2356 for ind in range(nChan):
2370 2357 data_param[ind,:,:] = self.__calculateMoments( spectra[ind,:,:] , absc , noise[ind] )
2371 #print data_param[:,3]
2372 2358 spwd = data_param[:,3]
2373 #print spwd.shape
2374 2359 # SPECB_JULIA,header=anal_header,jspectra=spectra,vel=velocities,hei=heights, num_aver=1, mode_fit=0,smoothing=smoothing,jvelr=velr,jspwd=spwd,jsnr=snr,jnoise=noise,jstdvnoise=stdvnoise
2375 #spwd1=[ 1.65607, 1.43416, 0.500373, 0.208361, 0.000000, 26.7767, 22.5936, 26.7530, 20.6962, 29.1098, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 28.0300, 27.0511, 27.8810, 26.3126, 27.8445, 24.6181, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000]
2376 #spwd=numpy.array([spwd1,spwd1,spwd1,spwd1])
2377 #print spwd.shape, heights.shape,coh_aver.shape
2378 2360 # para obtener spwd
2379 2361 for ic in range(nPairs):
2380 2362 pair = crosspairs[ic]
@@ -2399,23 +2381,15 class SpectralFitting(Operation):
2399 2381
2400 2382 return clean_coh_spectra, clean_coh_cspectra, clean_coh_aver
2401 2383
2402 #def CleanRayleigh(self,dataOut,spectra,cspectra,out_spectra,out_cspectra,sat_spectra,sat_cspectra,crosspairs,heights, channels, nProf,nHei,nChan,nPairs,nIncohInt,nBlocks):
2403 2384 def CleanRayleigh(self,dataOut,spectra,cspectra,save_drifts):
2404 #import matplotlib.pyplot as plt
2405 #for k in range(149):
2406
2407 # self.bloque0[:,:,:,k] = spectra[:,:,0:nHei]
2408 # self.bloques[:,:,:,k] = cspectra[:,:,0:nHei]
2409 #if self.i==nBlocks:
2410 # self.i==0
2411 rfunc = cspectra.copy() #self.bloques
2385
2386 rfunc = cspectra.copy()
2412 2387 n_funct = len(rfunc[0,:,0,0])
2413 val_spc = spectra*0.0 #self.bloque0*0.0
2414 val_cspc = cspectra*0.0 #self.bloques*0.0
2415 in_sat_spectra = spectra.copy() #self.bloque0
2416 in_sat_cspectra = cspectra.copy() #self.bloques
2388 val_spc = spectra*0.0
2389 val_cspc = cspectra*0.0
2390 in_sat_spectra = spectra.copy()
2391 in_sat_cspectra = cspectra.copy()
2417 2392
2418 #print( rfunc.shape)
2419 2393 min_hei = 200
2420 2394 nProf = dataOut.nProfiles
2421 2395 heights = dataOut.heightList
@@ -2426,22 +2400,19 class SpectralFitting(Operation):
2426 2400 nPairs = len(crosspairs)
2427 2401 hval=(heights >= min_hei).nonzero()
2428 2402 ih=hval[0]
2429 #print numpy.absolute(rfunc[:,0,0,14])
2430 2403 for ih in range(hval[0][0],nHei):
2431 2404 for ifreq in range(nProf):
2432 2405 for ii in range(n_funct):
2433 2406
2434 2407 func2clean = 10*numpy.log10(numpy.absolute(rfunc[:,ii,ifreq,ih]))
2435 #print numpy.amin(func2clean)
2436 2408 val = (numpy.isfinite(func2clean)==True).nonzero()
2437 2409 if len(val)>0:
2438 2410 min_val = numpy.around(numpy.amin(func2clean)-2) #> (-40)
2439 2411 if min_val <= -40 : min_val = -40
2440 2412 max_val = numpy.around(numpy.amax(func2clean)+2) #< 200
2441 2413 if max_val >= 200 : max_val = 200
2442 #print min_val, max_val
2443 2414 step = 1
2444 #Getting bins and the histogram
2415 #Getting bins and the histogram
2445 2416 x_dist = min_val + numpy.arange(1 + ((max_val-(min_val))/step))*step
2446 2417 y_dist,binstep = numpy.histogram(func2clean,bins=range(int(min_val),int(max_val+2),step))
2447 2418 mean = numpy.sum(x_dist * y_dist) / numpy.sum(y_dist)
@@ -2454,9 +2425,6 class SpectralFitting(Operation):
2454 2425 except:
2455 2426 mode = mean
2456 2427 stdv = sigma
2457 # 7.84616 53.9307 3.61863
2458 #stdv = 3.61863 # 2.99089
2459 #mode = 53.9307 #7.79008
2460 2428
2461 2429 #Removing echoes greater than mode + 3*stdv
2462 2430 factor_stdv = 2.5
@@ -2467,38 +2435,14 class SpectralFitting(Operation):
2467 2435 cross_pairs = crosspairs[ii]
2468 2436 #Getting coherent echoes which are removed.
2469 2437 if len(novall[0]) > 0:
2470 #val_spc[(0,1),novall[a],ih] = 1
2471 #val_spc[,(2,3),novall[a],ih] = 1
2472 2438 val_spc[novall[0],cross_pairs[0],ifreq,ih] = 1
2473 2439 val_spc[novall[0],cross_pairs[1],ifreq,ih] = 1
2474 2440 val_cspc[novall[0],ii,ifreq,ih] = 1
2475 #print("OUT NOVALL 1")
2476 2441 #Removing coherent from ISR data
2477 2442 spectra[noval,cross_pairs[0],ifreq,ih] = numpy.nan
2478 2443 spectra[noval,cross_pairs[1],ifreq,ih] = numpy.nan
2479 2444 cspectra[noval,ii,ifreq,ih] = numpy.nan
2480 #no sale es para savedrifts >2
2481 ''' channels = channels
2482 cross_pairs = cross_pairs
2483 #print("OUT NOVALL 2")
2484
2485 vcross0 = (cross_pairs[0] == channels[ii]).nonzero()
2486 vcross1 = (cross_pairs[1] == channels[ii]).nonzero()
2487 vcross = numpy.concatenate((vcross0,vcross1),axis=None)
2488 #print('vcros =', vcross)
2489
2490 #Getting coherent echoes which are removed.
2491 if len(novall) > 0:
2492 #val_spc[novall,ii,ifreq,ih] = 1
2493 val_spc[ii,ifreq,ih,novall] = 1
2494 if len(vcross) > 0:
2495 val_cspc[vcross,ifreq,ih,novall] = 1
2496
2497 #Removing coherent from ISR data.
2498 self.bloque0[ii,ifreq,ih,noval] = numpy.nan
2499 if len(vcross) > 0:
2500 self.bloques[vcross,ifreq,ih,noval] = numpy.nan
2501 '''
2445
2502 2446 #Getting average of the spectra and cross-spectra from incoherent echoes.
2503 2447 out_spectra = numpy.zeros([nChan,nProf,nHei], dtype=float) #+numpy.nan
2504 2448 out_cspectra = numpy.zeros([nPairs,nProf,nHei], dtype=complex) #+numpy.nan
@@ -2506,19 +2450,15 class SpectralFitting(Operation):
2506 2450 for ifreq in range(nProf):
2507 2451 for ich in range(nChan):
2508 2452 tmp = spectra[:,ich,ifreq,ih]
2509 valid = (numpy.isfinite(tmp[:])==True).nonzero()
2510 #print('TMP',tmp)
2453 valid = (numpy.isfinite(tmp[:])==True).nonzero()
2511 2454 if len(valid[0]) >0 :
2512 2455 out_spectra[ich,ifreq,ih] = numpy.nansum(tmp)/len(valid[0])
2513 #for icr in range(nPairs):
2514 2456 for icr in range(nPairs):
2515 2457 tmp = numpy.squeeze(cspectra[:,icr,ifreq,ih])
2516 2458 valid = (numpy.isfinite(tmp)==True).nonzero()
2517 2459 if len(valid[0]) > 0:
2518 2460 out_cspectra[icr,ifreq,ih] = numpy.nansum(tmp)/len(valid[0])
2519 # print('##########################################################')
2520 2461 #Removing fake coherent echoes (at least 4 points around the point)
2521
2522 2462 val_spectra = numpy.sum(val_spc,0)
2523 2463 val_cspectra = numpy.sum(val_cspc,0)
2524 2464
@@ -2564,39 +2504,26 class SpectralFitting(Operation):
2564 2504 valid = (numpy.isfinite(tmp)).nonzero()
2565 2505 if len(valid[0]) > 0:
2566 2506 sat_cspectra[icr,ifreq,ih] = numpy.nansum(tmp)/len(valid[0])
2567 #self.__dataReady= True
2568 #sat_spectra, sat_cspectra= sat_spectra, sat_cspectra
2569 #if not self.__dataReady:
2570 #return None, None
2571 2507 return out_spectra, out_cspectra,sat_spectra,sat_cspectra
2572 2508 def REM_ISOLATED_POINTS(self,array,rth):
2573 2509 if rth == None : rth = 4
2574
2575 2510 num_prof = len(array[0,:,0])
2576 2511 num_hei = len(array[0,0,:])
2577 2512 n2d = len(array[:,0,0])
2578 2513
2579 2514 for ii in range(n2d) :
2580 #print ii,n2d
2581 2515 tmp = array[ii,:,:]
2582 #indxs = WHERE(FINITE(tmp) AND tmp GT 0,cindxs)
2583 2516 tmp = numpy.reshape(tmp,num_prof*num_hei)
2584 2517 indxs1 = (numpy.isfinite(tmp)==True).nonzero()
2585 indxs2 = (tmp > 0).nonzero()
2586
2518 indxs2 = (tmp > 0).nonzero()
2587 2519 indxs1 = (indxs1[0])
2588 2520 indxs2 = indxs2[0]
2589 #indxs1 = numpy.array(indxs1[0])
2590 #indxs2 = numpy.array(indxs2[0])
2591 2521 indxs = None
2592 #print indxs1 , indxs2
2593 2522 for iv in range(len(indxs2)):
2594 2523 indv = numpy.array((indxs1 == indxs2[iv]).nonzero())
2595 #print len(indxs2), indv
2596 2524 if len(indv[0]) > 0 :
2597 2525 indxs = numpy.concatenate((indxs,indxs2[iv]), axis=None)
2598 2526 indxs = indxs[1:]
2599 #print indxs, len(indxs)
2600 2527 if len(indxs) < 4 :
2601 2528 array[ii,:,:] = 0.
2602 2529 return
@@ -2604,41 +2531,30 class SpectralFitting(Operation):
2604 2531 xpos = numpy.mod(indxs ,num_hei)
2605 2532 ypos = (indxs / num_hei)
2606 2533 sx = numpy.argsort(xpos) # Ordering respect to "x" (time)
2607 #print sx
2608 2534 xpos = xpos[sx]
2609 2535 ypos = ypos[sx]
2610 # *********************************** Cleaning isolated points **********************************
2536 # *********************************** Cleaning isolated points **********************************
2611 2537 ic = 0
2612 2538 while True :
2613 2539 r = numpy.sqrt(list(numpy.power((xpos[ic]-xpos),2)+ numpy.power((ypos[ic]-ypos),2)))
2614 #no_coh = WHERE(FINITE(r) AND (r LE rth),cno_coh)
2615 #plt.plot(r)
2616 #plt.show()
2617 2540 no_coh1 = (numpy.isfinite(r)==True).nonzero()
2618 2541 no_coh2 = (r <= rth).nonzero()
2619 #print r, no_coh1, no_coh2
2620 2542 no_coh1 = numpy.array(no_coh1[0])
2621 2543 no_coh2 = numpy.array(no_coh2[0])
2622 2544 no_coh = None
2623 #print valid1 , valid2
2624 2545 for iv in range(len(no_coh2)):
2625 2546 indv = numpy.array((no_coh1 == no_coh2[iv]).nonzero())
2626 2547 if len(indv[0]) > 0 :
2627 2548 no_coh = numpy.concatenate((no_coh,no_coh2[iv]), axis=None)
2628 2549 no_coh = no_coh[1:]
2629 #print len(no_coh), no_coh
2630 2550 if len(no_coh) < 4 :
2631 #print xpos[ic], ypos[ic], ic
2632 2551 xpos[ic] = numpy.nan
2633 2552 ypos[ic] = numpy.nan
2634 2553
2635 2554 ic = ic + 1
2636 2555 if (ic == len(indxs)) :
2637 2556 break
2638 #print( xpos, ypos)
2639
2640 2557 indxs = (numpy.isfinite(list(xpos))==True).nonzero()
2641 #print indxs[0]
2642 2558 if len(indxs[0]) < 4 :
2643 2559 array[ii,:,:] = 0.
2644 2560 return
@@ -2652,20 +2568,12 class SpectralFitting(Operation):
2652 2568
2653 2569 tmp[list(xpos + (ypos*num_hei))] = junk[list(xpos + (ypos*num_hei))]
2654 2570 array[ii,:,:] = numpy.reshape(tmp,(num_prof,num_hei))
2655
2656 #print array.shape
2657 #tmp = numpy.reshape(tmp,(num_prof,num_hei))
2658 #print tmp.shape
2659
2660 2571 return array
2572
2661 2573 def moments(self,doppler,yarray,npoints):
2662 2574 ytemp = yarray
2663 #val = WHERE(ytemp GT 0,cval)
2664 #if cval == 0 : val = range(npoints-1)
2665 2575 val = (ytemp > 0).nonzero()
2666 2576 val = val[0]
2667 #print('hvalid:',hvalid)
2668 #print('valid', valid)
2669 2577 if len(val) == 0 : val = range(npoints-1)
2670 2578
2671 2579 ynew = 0.5*(ytemp[val[0]]+ytemp[val[len(val)-1]])
@@ -2684,111 +2592,170 class SpectralFitting(Operation):
2684 2592
2685 2593
2686 2594
2595 def run(self, dataOut, getSNR = True, path=None, file=None, groupList=None, filec=None,coh_th=None, hei_th=None,taver=None,proc=None,nhei=None,nprofs=None,ipp=None,channelList=None):
2596 if not numpy.any(proc):
2597 nChannels = dataOut.nChannels
2598 nHeights= dataOut.heightList.size
2599 nProf = dataOut.nProfiles
2600 if numpy.any(taver): taver=int(taver)
2601 else : taver = 5
2602 tini=time.localtime(dataOut.utctime)
2603 if (tini.tm_min % taver) == 0 and (tini.tm_sec < 5 and self.fint==0):
2604 self.index = 0
2605 jspc = self.buffer
2606 jcspc = self.buffer2
2607 jnoise = self.buffer3
2608 self.buffer = dataOut.data_spc
2609 self.buffer2 = dataOut.data_cspc
2610 self.buffer3 = dataOut.noise
2611 self.fint = 1
2612 if numpy.any(jspc) :
2613 jspc= numpy.reshape(jspc,(int(len(jspc)/nChannels),nChannels,nProf,nHeights))
2614 jcspc= numpy.reshape(jcspc,(int(len(jcspc)/int(nChannels/2)),int(nChannels/2),nProf,nHeights))
2615 jnoise= numpy.reshape(jnoise,(int(len(jnoise)/nChannels),nChannels))
2616 else:
2617 dataOut.flagNoData = True
2618 return dataOut
2619 else :
2620 if (tini.tm_min % taver) == 0 : self.fint = 1
2621 else : self.fint = 0
2622 self.index += 1
2623 if numpy.any(self.buffer):
2624 self.buffer = numpy.concatenate((self.buffer,dataOut.data_spc), axis=0)
2625 self.buffer2 = numpy.concatenate((self.buffer2,dataOut.data_cspc), axis=0)
2626 self.buffer3 = numpy.concatenate((self.buffer3,dataOut.noise), axis=0)
2627 else:
2628 self.buffer = dataOut.data_spc
2629 self.buffer2 = dataOut.data_cspc
2630 self.buffer3 = dataOut.noise
2631 dataOut.flagNoData = True
2632 return dataOut
2633 if path != None:
2634 sys.path.append(path)
2635 self.library = importlib.import_module(file)
2636 if filec != None:
2637 self.weightf = importlib.import_module(filec)
2638
2639 #To be inserted as a parameter
2640 groupArray = numpy.array(groupList)
2641 #groupArray = numpy.array([[0,1],[2,3]])
2642 dataOut.groupList = groupArray
2643 nGroups = groupArray.shape[0]
2644 nChannels = dataOut.nChannels
2645 nHeights = dataOut.heightList.size
2646
2647 #Parameters Array
2648 dataOut.data_param = None
2649 dataOut.data_paramC = None
2650 dataOut.clean_num_aver = None
2651 dataOut.coh_num_aver = None
2652 dataOut.tmp_spectra_i = None
2653 dataOut.tmp_cspectra_i = None
2654 dataOut.tmp_spectra_c = None
2655 dataOut.tmp_cspectra_c = None
2656 dataOut.index = None
2657
2658 #Set constants
2659 constants = self.library.setConstants(dataOut)
2660 dataOut.constants = constants
2661 M = dataOut.normFactor
2662 N = dataOut.nFFTPoints
2663 ippSeconds = dataOut.ippSeconds
2664 K = dataOut.nIncohInt
2665 pairsArray = numpy.array(dataOut.pairsList)
2666 snrth= 20
2667 spectra = dataOut.data_spc
2668 cspectra = dataOut.data_cspc
2669 nProf = dataOut.nProfiles
2670 heights = dataOut.heightList
2671 nHei = len(heights)
2672 channels = dataOut.channelList
2673 nChan = len(channels)
2674 nIncohInt = dataOut.nIncohInt
2675 crosspairs = dataOut.groupList
2676 noise = dataOut.noise
2677 jnoise = jnoise/N
2678 noise = numpy.nansum(jnoise,axis=0)#/len(jnoise)
2679 power = numpy.sum(spectra, axis=1)
2680 nPairs = len(crosspairs)
2681 absc = dataOut.abscissaList[:-1]
2687 2682
2688
2689
2690 def run(self, dataOut, getSNR = True, path=None, file=None, groupList=None, filec=None,coh_th=None, hei_th=None):
2691 nChannels = dataOut.nChannels
2692 nHeights= dataOut.heightList.size
2693 nProf = dataOut.nProfiles
2694 tini=time.localtime(dataOut.utctime)
2695 if (tini.tm_min % 5) == 0 and (tini.tm_sec < 5 and self.fint==0):
2696 self.index = 0
2697 jspc = self.buffer
2698 jcspc = self.buffer2
2699 jnoise = self.buffer3
2700 self.buffer = dataOut.data_spc
2701 self.buffer2 = dataOut.data_cspc
2702 self.buffer3 = dataOut.noise
2703 self.fint = 1
2704 if numpy.any(jspc) :
2705 jspc= numpy.reshape(jspc,(int(len(jspc)/4),nChannels,nProf,nHeights))
2706 jcspc= numpy.reshape(jcspc,(int(len(jcspc)/2),2,nProf,nHeights))
2707 jnoise= numpy.reshape(jnoise,(int(len(jnoise)/4),nChannels))
2708 else:
2709 dataOut.flagNoData = True
2710 return dataOut
2711 else :
2712 if (tini.tm_min % 5) == 0 : self.fint = 1
2713 else : self.fint = 0
2714 self.index += 1
2715 if numpy.any(self.buffer):
2716 self.buffer = numpy.concatenate((self.buffer,dataOut.data_spc), axis=0)
2717 self.buffer2 = numpy.concatenate((self.buffer2,dataOut.data_cspc), axis=0)
2718 self.buffer3 = numpy.concatenate((self.buffer3,dataOut.noise), axis=0)
2719 else:
2720 self.buffer = dataOut.data_spc
2721 self.buffer2 = dataOut.data_cspc
2722 self.buffer3 = dataOut.noise
2723 dataOut.flagNoData = True
2724 return dataOut
2725 if path != None:
2726 sys.path.append(path)
2727 self.library = importlib.import_module(file)
2728
2729 #To be inserted as a parameter
2683 if not self.isConfig:
2684 self.isConfig = True
2685
2686 index = tini.tm_hour*12+tini.tm_min/taver
2687 dataOut.index= index
2688 jspc = jspc/N/N
2689 jcspc = jcspc/N/N
2690 tmp_spectra,tmp_cspectra,sat_spectra,sat_cspectra = self.CleanRayleigh(dataOut,jspc,jcspc,2)
2691 jspectra = tmp_spectra*len(jspc[:,0,0,0])
2692 jcspectra = tmp_cspectra*len(jspc[:,0,0,0])
2693 my_incoh_spectra ,my_incoh_cspectra,my_incoh_aver,my_coh_aver, incoh_spectra, coh_spectra, incoh_cspectra, coh_cspectra, incoh_aver, coh_aver = self.__DiffCoherent(jspectra, jcspectra, dataOut, noise, snrth,coh_th, hei_th)
2694 clean_coh_spectra, clean_coh_cspectra, clean_coh_aver = self.__CleanCoherent(snrth, coh_spectra, coh_cspectra, coh_aver, dataOut, noise,1,index)
2695 dataOut.data_spc = incoh_spectra
2696 dataOut.data_cspc = incoh_cspectra
2697 clean_num_aver = incoh_aver*len(jspc[:,0,0,0])
2698 coh_num_aver = clean_coh_aver*len(jspc[:,0,0,0])
2699 dataOut.clean_num_aver = clean_num_aver
2700 dataOut.coh_num_aver = coh_num_aver
2701 dataOut.tmp_spectra_i = incoh_spectra
2702 dataOut.tmp_cspectra_i = incoh_cspectra
2703 dataOut.tmp_spectra_c = clean_coh_spectra
2704 dataOut.tmp_cspectra_c = clean_coh_cspectra
2705 #List of possible combinations
2706 listComb = itertools.combinations(numpy.arange(groupArray.shape[1]),2)
2707 indCross = numpy.zeros(len(list(listComb)), dtype = 'int')
2708 if getSNR:
2709 listChannels = groupArray.reshape((groupArray.size))
2710 listChannels.sort()
2711 dataOut.data_SNR = self.__getSNR(dataOut.data_spc[listChannels,:,:], noise[listChannels])
2712 else:
2713 clean_num_aver = dataOut.clean_num_aver
2714 coh_num_aver = dataOut.coh_num_aver
2715 dataOut.data_spc = dataOut.tmp_spectra_i
2716 dataOut.data_cspc = dataOut.tmp_cspectra_i
2717 clean_coh_spectra = dataOut.tmp_spectra_c
2718 clean_coh_cspectra = dataOut.tmp_cspectra_c
2719 jspectra = dataOut.data_spc+clean_coh_spectra
2720 nHeights = len(dataOut.heightList) # nhei
2721 nProf = int(dataOut.nProfiles)
2722 dataOut.nProfiles = nProf
2723 dataOut.data_param = None
2724 dataOut.data_paramC = None
2725 dataOut.code = numpy.array([[-1.,-1.,1.],[1.,1.,-1.]])
2726 #M=600
2727 #N=200
2728 dataOut.flagDecodeData=True
2729 M = int(dataOut.normFactor)
2730 N = int(dataOut.nFFTPoints)
2731 dataOut.nFFTPoints = N
2732 dataOut.nIncohInt= int(dataOut.nIncohInt)
2733 dataOut.nProfiles = int(dataOut.nProfiles)
2734 dataOut.nCohInt = int(dataOut.nCohInt)
2735 print('sale',dataOut.nProfiles,dataOut.nHeights)
2736 #dataOut.nFFTPoints=nprofs
2737 #dataOut.normFactor = nprofs
2738 dataOut.channelList = channelList
2739 #dataOut.ippFactor=1
2740 #ipp = ipp/150*1.e-3
2741 vmax = (300000000/49920000.0/2) / (dataOut.ippSeconds)
2742 #dataOut.ippSeconds=ipp
2743 absc = vmax*( numpy.arange(nProf,dtype='float')-nProf/2.)/nProf
2744 print('sale 2',dataOut.ippSeconds,M,N)
2745 print('Empieza procesamiento offline')
2746 if path != None:
2747 sys.path.append(path)
2748 self.library = importlib.import_module(file)
2749 constants = self.library.setConstants(dataOut)
2750 constants['M'] = M
2751 dataOut.constants = constants
2752
2730 2753 groupArray = numpy.array(groupList)
2731 #groupArray = numpy.array([[0,1],[2,3]])
2732 2754 dataOut.groupList = groupArray
2733
2734 2755 nGroups = groupArray.shape[0]
2735 nChannels = dataOut.nChannels
2736 nHeights = dataOut.heightList.size
2737
2738 #Parameters Array
2739 dataOut.data_param = None
2740 dataOut.data_paramC = None
2741
2742 #Set constants
2743 constants = self.library.setConstants(dataOut)
2744 dataOut.constants = constants
2745 M = dataOut.normFactor
2746 N = dataOut.nFFTPoints
2747 ippSeconds = dataOut.ippSeconds
2748 K = dataOut.nIncohInt
2749 pairsArray = numpy.array(dataOut.pairsList)
2750
2751 snrth= 20
2752 spectra = dataOut.data_spc
2753 cspectra = dataOut.data_cspc
2754 nProf = dataOut.nProfiles
2755 heights = dataOut.heightList
2756 nHei = len(heights)
2757 channels = dataOut.channelList
2758 nChan = len(channels)
2759 nIncohInt = dataOut.nIncohInt
2760 crosspairs = dataOut.groupList
2761 noise = dataOut.noise
2762 jnoise = jnoise/N
2763 noise = numpy.nansum(jnoise,axis=0)#/len(jnoise)
2764 power = numpy.sum(spectra, axis=1)
2765 nPairs = len(crosspairs)
2766 absc = dataOut.abscissaList[:-1]
2767
2768 if not self.isConfig:
2769 self.isConfig = True
2770
2771 index = tini.tm_hour*12+tini.tm_min/5
2772 jspc = jspc/N/N
2773 jcspc = jcspc/N/N
2774 tmp_spectra,tmp_cspectra,sat_spectra,sat_cspectra = self.CleanRayleigh(dataOut,jspc,jcspc,2)
2775 jspectra = tmp_spectra*len(jspc[:,0,0,0])
2776 jcspectra = tmp_cspectra*len(jspc[:,0,0,0])
2777 my_incoh_spectra ,my_incoh_cspectra,my_incoh_aver,my_coh_aver, incoh_spectra, coh_spectra, incoh_cspectra, coh_cspectra, incoh_aver, coh_aver = self.__DiffCoherent(jspectra, jcspectra, dataOut, noise, snrth, coh_th, hei_th)
2778 clean_coh_spectra, clean_coh_cspectra, clean_coh_aver = self.__CleanCoherent(snrth, coh_spectra, coh_cspectra, coh_aver, dataOut, noise,1,index)
2779 dataOut.data_spc = incoh_spectra
2780 dataOut.data_cspc = incoh_cspectra
2781
2782 clean_num_aver = incoh_aver*len(jspc[:,0,0,0])
2783 coh_num_aver = clean_coh_aver*len(jspc[:,0,0,0])
2784 2756 #List of possible combinations
2785 2757 listComb = itertools.combinations(numpy.arange(groupArray.shape[1]),2)
2786 2758 indCross = numpy.zeros(len(list(listComb)), dtype = 'int')
2787
2788 if getSNR:
2789 listChannels = groupArray.reshape((groupArray.size))
2790 listChannels.sort()
2791 dataOut.data_SNR = self.__getSNR(dataOut.data_spc[listChannels,:,:], noise[listChannels])
2792 2759 if dataOut.data_paramC is None:
2793 2760 dataOut.data_paramC = numpy.zeros((nGroups*4, nHeights,2))*numpy.nan
2794 2761 for i in range(nGroups):
@@ -2807,17 +2774,18 class SpectralFitting(Operation):
2807 2774 dataCross = dataCross**2
2808 2775 nhei = nHeights
2809 2776 poweri = numpy.sum(dataOut.data_spc[:,1:nProf-0,:],axis=1)/clean_num_aver[:,:]
2810 if i == 0 : my_noises = numpy.zeros(4,dtype=float) #FLTARR(4)
2777 if i == 0 : my_noises = numpy.zeros(4,dtype=float)
2811 2778 n0i = numpy.nanmin(poweri[0+i*2,0:nhei-0])/(nProf-1)
2812 2779 n1i = numpy.nanmin(poweri[1+i*2,0:nhei-0])/(nProf-1)
2813 2780 n0 = n0i
2814 2781 n1= n1i
2815 2782 my_noises[2*i+0] = n0
2816 2783 my_noises[2*i+1] = n1
2817 snrth = -15.0 #-16 -10
2784 snrth = -15.0 # -4 -16 -25
2818 2785 snrth = 10**(snrth/10.0)
2819 2786 jvelr = numpy.zeros(nHeights, dtype = 'float')
2820 2787 hvalid = [0]
2788 coh2 = abs(dataOut.data_cspc[i,1:nProf,:])**2/(dataOut.data_spc[0+i*2,1:nProf-0,:]*dataOut.data_spc[1+i*2,1:nProf-0,:])
2821 2789 for h in range(nHeights):
2822 2790 smooth = clean_num_aver[i+1,h]
2823 2791 signalpn0 = (dataOut.data_spc[i*2,1:(nProf-0),h])/smooth
@@ -2826,15 +2794,23 class SpectralFitting(Operation):
2826 2794 signal1 = signalpn1-n1
2827 2795 snr0 = numpy.sum(signal0/n0)/(nProf-1)
2828 2796 snr1 = numpy.sum(signal1/n1)/(nProf-1)
2829 maxp0 = numpy.argmax(signal0)
2830 maxp1 = numpy.argmax(signal1)
2831 jvelr[h] = (absc[maxp0]+absc[maxp1])/2.
2797 gamma = coh2[:,h]
2798 indxs = (numpy.isfinite(list(gamma))==True).nonzero()
2799 if len(indxs) >0:
2800 if numpy.nanmean(gamma) > 0.07:
2801 maxp0 = numpy.argmax(signal0*gamma)
2802 maxp1 = numpy.argmax(signal1*gamma)
2803 #print('usa gamma',numpy.nanmean(gamma))
2804 else:
2805 maxp0 = numpy.argmax(signal0)
2806 maxp1 = numpy.argmax(signal1)
2807 jvelr[h] = (absc[maxp0]+absc[maxp1])/2.
2808 else: jvelr[h] = absc[0]
2832 2809 if snr0 > 0.1 and snr1 > 0.1: hvalid = numpy.concatenate((hvalid,h), axis=None)
2833 2810 #print(maxp0,absc[maxp0],snr0,jvelr[h])
2834 2811
2835 2812 if len(hvalid)> 1: fd0 = numpy.median(jvelr[hvalid[1:]])*-1
2836 2813 else: fd0 = numpy.nan
2837
2838 2814 for h in range(nHeights):
2839 2815 d = data[:,h]
2840 2816 smooth = clean_num_aver[i+1,h] #dataOut.data_spc[:,1:nProf-0,:]
@@ -2859,11 +2835,8 class SpectralFitting(Operation):
2859 2835 D[x*N:(x+1)*N,y*N:(y+1)*N] = D12
2860 2836 D[y*N:(y+1)*N,x*N:(x+1)*N] = D12
2861 2837 ind += 1
2862 #diagD = numpy.zeros(256)
2863 #if h == 17 :
2864 # for ii in range(256): diagD[ii] = D[ii,ii]
2865 #Dinv=numpy.linalg.inv(D)
2866 #L=numpy.linalg.cholesky(Dinv)
2838 diagD = numpy.zeros(256)
2839
2867 2840 try:
2868 2841 Dinv=numpy.linalg.inv(D)
2869 2842 L=numpy.linalg.cholesky(Dinv)
@@ -2873,14 +2846,18 class SpectralFitting(Operation):
2873 2846 LT=L.T
2874 2847
2875 2848 dp = numpy.dot(LT,d)
2876
2877 2849 #Initial values
2878 2850 data_spc = dataOut.data_spc[coord,:,h]
2851 w = data_spc/data_spc
2852 if filec != None:
2853 w = self.weightf.weightfit(w,tini.tm_year,tini.tm_yday,index,h,i)
2879 2854 if (h>6)and(error1[3]<25):
2880 2855 p0 = dataOut.data_param[i,:,h-1]
2881 2856 else:
2882 p0 = numpy.array(self.library.initialValuesFunction(data_spc, 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)
2883 2858 p0[3] = fd0
2859 if filec != None:
2860 p0 = self.weightf.Vrfit(p0,tini.tm_year,tini.tm_yday,index,h,i)
2884 2861 try:
2885 2862 #Least Squares
2886 2863 minp,covp,infodict,mesg,ier = optimize.leastsq(self.__residFunction,p0,args=(dp,LT,constants),full_output=True)
@@ -2899,7 +2876,7 class SpectralFitting(Operation):
2899 2876 p0 = numpy.array(self.library.initialValuesFunction(data_spc, constants))
2900 2877 minp = p0*numpy.nan
2901 2878 error0 = numpy.nan
2902 error1 = p0*numpy.nan
2879 error1 = p0*numpy.nan
2903 2880 if dataOut.data_param is None:
2904 2881 dataOut.data_param = numpy.zeros((nGroups, p0.size, nHeights))*numpy.nan
2905 2882 dataOut.data_error = numpy.zeros((nGroups, p0.size + 1, nHeights))*numpy.nan
@@ -2924,11 +2901,11 class SpectralFitting(Operation):
2924 2901 dataOut.data_paramC[0+4*i,ht,0] = numpy.sum((signalpn0/val0_npoints))/n0
2925 2902 dataOut.data_paramC[1+4*i,ht,0] = numpy.sum((signalpn1/val1_npoints))/n1
2926 2903
2927 signal0 = (signalpn0-n0) # > 0
2904 signal0 = (signalpn0-n0)
2928 2905 vali = (signal0 < 0).nonzero()
2929 2906 vali = vali[0]
2930 2907 if len(vali) > 0 : signal0[vali] = 0
2931 signal1 = (signalpn1-n1) #> 0
2908 signal1 = (signalpn1-n1)
2932 2909 vali = (signal1 < 0).nonzero()
2933 2910 vali = vali[0]
2934 2911 if len(vali) > 0 : signal1[vali] = 0
@@ -2947,6 +2924,7 class SpectralFitting(Operation):
2947 2924
2948 2925 dataOut.data_spc = jspectra
2949 2926 dataOut.spc_noise = my_noises*nProf*M
2927 if numpy.any(proc): dataOut.spc_noise = my_noises*nProf*M
2950 2928 if getSNR:
2951 2929 listChannels = groupArray.reshape((groupArray.size))
2952 2930 listChannels.sort()
@@ -3468,8 +3446,8 class WindProfiler(Operation):
3468 3446
3469 3447 def run(self, dataOut, technique, nHours=1, hmin=70, hmax=110, **kwargs):
3470 3448
3471 param = dataOut.data_param
3472 if dataOut.abscissaList != None:
3449 param = dataOut.moments
3450 if numpy.any(dataOut.abscissaList):
3473 3451 absc = dataOut.abscissaList[:-1]
3474 3452 # noise = dataOut.noise
3475 3453 heightList = dataOut.heightList
@@ -3610,7 +3588,7 class WindProfiler(Operation):
3610 3588 dataOut.flagNoData = False
3611 3589 self.__buffer = None
3612 3590
3613 return
3591 return dataOut
3614 3592
3615 3593 class EWDriftsEstimation(Operation):
3616 3594
@@ -3638,18 +3616,13 class EWDriftsEstimation(Operation):
3638 3616 y1=y1[vali]
3639 3617 x = x[vali]
3640 3618 f1 = interpolate.interp1d(x,y1,kind = 'cubic',bounds_error=False)
3641
3642 #heiRang1 = x*math.cos(phi[maxid])
3643 3619 x1 = heiRang1
3644 3620 y11 = f1(x1)
3645
3646 3621 y2 = SNR[i,:]
3647 #print 'snr ', y2
3648 3622 x = heiRang*math.cos(phi[i])
3649 3623 vali= (y2 != -1).nonzero()
3650 3624 y2 = y2[vali]
3651 3625 x = x[vali]
3652 #print 'snr ',y2
3653 3626 f2 = interpolate.interp1d(x,y2,kind = 'cubic',bounds_error=False)
3654 3627 y21 = f2(x1)
3655 3628
@@ -3663,12 +3636,10 class EWDriftsEstimation(Operation):
3663 3636 heiRang = dataOut.heightList
3664 3637 velRadial = dataOut.data_param[:,3,:]
3665 3638 velRadialm = dataOut.data_param[:,2:4,:]*-1
3666
3667 3639 rbufc=dataOut.data_paramC[:,:,0]
3668 3640 ebufc=dataOut.data_paramC[:,:,1]
3669 3641 SNR = dataOut.data_snr
3670 3642 velRerr = dataOut.data_error[:,4,:]
3671 #moments=numpy.vstack(([velRadialm[0,:]],[velRadialm[0,:]],[velRadialm[1,:]],[velRadialm[1,:]]))
3672 3643 channels = dataOut.channelList
3673 3644 nChan = len(channels)
3674 3645 my_nbeams = nChan/2
@@ -3683,7 +3654,6 class EWDriftsEstimation(Operation):
3683 3654 p_w1C = rbufc[1,:]
3684 3655 w_wC = rbufc[2,:]*-1 #*radial_sign(radial EQ 1)
3685 3656 t_wC = rbufc[3,:]
3686 #my_nbeams = 2
3687 3657 if my_nbeams == 1:
3688 3658 w = velRadial[0,:]
3689 3659 winds = velRadial.copy()
@@ -3714,7 +3684,6 class EWDriftsEstimation(Operation):
3714 3684 if len(bad) > 0 :
3715 3685 w_w[bad] = w_wC[bad]
3716 3686 w_w_err[bad]= numpy.nan
3717 # if my_nbeams == 2:
3718 3687 smooth_eC=ebufc[4,:]
3719 3688 p_e0C = rbufc[4,:]
3720 3689 p_e1C = rbufc[5,:]
@@ -3742,17 +3711,15 class EWDriftsEstimation(Operation):
3742 3711 dataOut.data_snr1 = numpy.reshape(snr1,(1,snr1.shape[0]))
3743 3712 dataOut.utctimeInit = dataOut.utctime
3744 3713 dataOut.outputInterval = dataOut.timeInterval
3745
3714
3746 3715 hei_aver0 = 218
3747 3716 jrange = 450 #900 para HA drifts
3748 deltah = 15 # para EWDrifts 15.0 MP 15 #dataOut.spacing(0)
3717 deltah = 15.0 #dataOut.spacing(0) 25 HAD
3749 3718 h0 = 0.0 #dataOut.first_height(0)
3750 3719 heights = dataOut.heightList
3751 3720 nhei = len(heights)
3752 3721
3753 3722 range1 = numpy.arange(nhei) * deltah + h0
3754
3755 #jhei = WHERE(range1 GE hei_aver0 , jcount)
3756 3723 jhei = (range1 >= hei_aver0).nonzero()
3757 3724 if len(jhei[0]) > 0 :
3758 3725 h0_index = jhei[0][0] # Initial height for getting averages 218km
@@ -3776,7 +3743,6 class EWDriftsEstimation(Operation):
3776 3743
3777 3744 wA = w[h0_index:h0_index+nhei_avg-0]
3778 3745 wA_err = w_err[h0_index:h0_index+nhei_avg-0]
3779
3780 3746 for i in range(5) :
3781 3747 vals = wA[i*h_avgs:(i+1)*h_avgs-0]
3782 3748 errs = wA_err[i*h_avgs:(i+1)*h_avgs-0]
@@ -3796,7 +3762,6 class EWDriftsEstimation(Operation):
3796 3762 wA = wA[6*h_avgs:nhei_avg-0]
3797 3763 wA_err=wA_err[6*h_avgs:nhei_avg-0]
3798 3764 if my_nbeams == 2 :
3799
3800 3765 uA = u[h0_index:h0_index+nhei_avg]
3801 3766 uA_err=u_err[h0_index:h0_index+nhei_avg]
3802 3767
@@ -3821,12 +3786,9 class EWDriftsEstimation(Operation):
3821 3786 tini=time.localtime(dataOut.utctime)
3822 3787 datefile= str(tini[0]).zfill(4)+str(tini[1]).zfill(2)+str(tini[2]).zfill(2)
3823 3788 nfile = fileDrifts+'/jro'+datefile+'drifts.txt'
3824
3825 3789 f1 = open(nfile,'a')
3826
3827 3790 datedriftavg=str(tini[0])+' '+str(tini[1])+' '+str(tini[2])+' '+str(tini[3])+' '+str(tini[4])
3828 3791 driftavgstr=str(dataOut.drifts_avg)
3829
3830 3792 numpy.savetxt(f1,numpy.column_stack([tini[0],tini[1],tini[2],tini[3],tini[4]]),fmt='%4i')
3831 3793 numpy.savetxt(f1,dataOut.drifts_avg,fmt='%10.2f')
3832 3794 f1.close()
@@ -5242,4 +5204,138 class IGRFModel(Operation):
5242 5204
5243 5205 mkfact_short_2020.mkfact(dataOut.year,dataOut.h,dataOut.bfm,dataOut.thb,dataOut.bki,dataOut.MAXNRANGENDT)
5244 5206
5245 return dataOut No newline at end of file
5207 return dataOut
5208
5209 class MergeProc(ProcessingUnit):
5210
5211 def __init__(self):
5212 ProcessingUnit.__init__(self)
5213
5214 def run(self, attr_data, attr_data_2 = None, attr_data_3 = None, attr_data_4 = None, attr_data_5 = None, mode=0):
5215
5216 self.dataOut = getattr(self, self.inputs[0])
5217 data_inputs = [getattr(self, attr) for attr in self.inputs]
5218 #print(self.inputs)
5219 #print(numpy.shape([getattr(data, attr_data) for data in data_inputs][1]))
5220 #exit(1)
5221 if mode==0:
5222 data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs])
5223 setattr(self.dataOut, attr_data, data)
5224
5225 if mode==1: #Hybrid
5226 #data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs],axis=1)
5227 #setattr(self.dataOut, attr_data, data)
5228 setattr(self.dataOut, 'dataLag_spc', [getattr(data, attr_data) for data in data_inputs][0])
5229 setattr(self.dataOut, 'dataLag_spc_LP', [getattr(data, attr_data) for data in data_inputs][1])
5230 setattr(self.dataOut, 'dataLag_cspc', [getattr(data, attr_data_2) for data in data_inputs][0])
5231 setattr(self.dataOut, 'dataLag_cspc_LP', [getattr(data, attr_data_2) for data in data_inputs][1])
5232 #setattr(self.dataOut, 'nIncohInt', [getattr(data, attr_data_3) for data in data_inputs][0])
5233 #setattr(self.dataOut, 'nIncohInt_LP', [getattr(data, attr_data_3) for data in data_inputs][1])
5234 '''
5235 print(self.dataOut.dataLag_spc_LP.shape)
5236 print(self.dataOut.dataLag_cspc_LP.shape)
5237 exit(1)
5238 '''
5239
5240 #self.dataOut.dataLag_spc_LP = numpy.transpose(self.dataOut.dataLag_spc_LP[0],(2,0,1))
5241 #self.dataOut.dataLag_cspc_LP = numpy.transpose(self.dataOut.dataLag_cspc_LP,(3,1,2,0))
5242 '''
5243 print("Merge")
5244 print(numpy.shape(self.dataOut.dataLag_spc))
5245 print(numpy.shape(self.dataOut.dataLag_spc_LP))
5246 print(numpy.shape(self.dataOut.dataLag_cspc))
5247 print(numpy.shape(self.dataOut.dataLag_cspc_LP))
5248 exit(1)
5249 '''
5250 #print(numpy.sum(self.dataOut.dataLag_spc_LP[2,:,164])/128)
5251 #print(numpy.sum(self.dataOut.dataLag_cspc_LP[0,:,30,1])/128)
5252 #exit(1)
5253 #print(self.dataOut.NDP)
5254 #print(self.dataOut.nNoiseProfiles)
5255
5256 #self.dataOut.nIncohInt_LP = 128
5257 self.dataOut.nProfiles_LP = 128#self.dataOut.nIncohInt_LP
5258 self.dataOut.nIncohInt_LP = self.dataOut.nIncohInt
5259 self.dataOut.NLAG = 16
5260 self.dataOut.NRANGE = 200
5261 self.dataOut.NSCAN = 128
5262 #print(numpy.shape(self.dataOut.data_spc))
5263
5264 #exit(1)
5265
5266 if mode==2: #HAE 2022
5267 data = numpy.sum([getattr(data, attr_data) for data in data_inputs],axis=0)
5268 setattr(self.dataOut, attr_data, data)
5269
5270 self.dataOut.nIncohInt *= 2
5271 #meta = self.dataOut.getFreqRange(1)/1000.
5272 self.dataOut.freqRange = self.dataOut.getFreqRange(1)/1000.
5273
5274 #exit(1)
5275
5276 if mode==4: #Hybrid LP-SSheightProfiles
5277 #data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs],axis=1)
5278 #setattr(self.dataOut, attr_data, data)
5279 setattr(self.dataOut, 'dataLag_spc', getattr(data_inputs[0], attr_data)) #DP
5280 setattr(self.dataOut, 'dataLag_cspc', getattr(data_inputs[0], attr_data_2)) #DP
5281 setattr(self.dataOut, 'dataLag_spc_LP', getattr(data_inputs[1], attr_data_3)) #LP
5282 #setattr(self.dataOut, 'dataLag_cspc_LP', getattr(data_inputs[1], attr_data_4)) #LP
5283 #setattr(self.dataOut, 'data_acf', getattr(data_inputs[1], attr_data_5)) #LP
5284 setattr(self.dataOut, 'data_acf', getattr(data_inputs[1], attr_data_5)) #LP
5285 #print("Merge data_acf: ",self.dataOut.data_acf.shape)
5286 #exit(1)
5287 #print(self.dataOut.data_spc_LP.shape)
5288 #print("Exit")
5289 #exit(1)
5290 #setattr(self.dataOut, 'dataLag_cspc', [getattr(data, attr_data_2) for data in data_inputs][0])
5291 #setattr(self.dataOut, 'dataLag_cspc_LP', [getattr(data, attr_data_2) for data in data_inputs][1])
5292 #setattr(self.dataOut, 'nIncohInt', [getattr(data, attr_data_3) for data in data_inputs][0])
5293 #setattr(self.dataOut, 'nIncohInt_LP', [getattr(data, attr_data_3) for data in data_inputs][1])
5294 '''
5295 print(self.dataOut.dataLag_spc_LP.shape)
5296 print(self.dataOut.dataLag_cspc_LP.shape)
5297 exit(1)
5298 '''
5299 '''
5300 print(self.dataOut.dataLag_spc_LP[0,:,100])
5301 print(self.dataOut.dataLag_spc_LP[1,:,100])
5302 exit(1)
5303 '''
5304 #self.dataOut.dataLag_spc_LP = numpy.transpose(self.dataOut.dataLag_spc_LP[0],(2,0,1))
5305 #self.dataOut.dataLag_cspc_LP = numpy.transpose(self.dataOut.dataLag_cspc_LP,(3,1,2,0))
5306 '''
5307 print("Merge")
5308 print(numpy.shape(self.dataOut.dataLag_spc))
5309 print(numpy.shape(self.dataOut.dataLag_spc_LP))
5310 print(numpy.shape(self.dataOut.dataLag_cspc))
5311 print(numpy.shape(self.dataOut.dataLag_cspc_LP))
5312 exit(1)
5313 '''
5314 #print(numpy.sum(self.dataOut.dataLag_spc_LP[2,:,164])/128)
5315 #print(numpy.sum(self.dataOut.dataLag_cspc_LP[0,:,30,1])/128)
5316 #exit(1)
5317 #print(self.dataOut.NDP)
5318 #print(self.dataOut.nNoiseProfiles)
5319
5320 #self.dataOut.nIncohInt_LP = 128
5321 #self.dataOut.nProfiles_LP = 128#self.dataOut.nIncohInt_LP
5322 self.dataOut.nProfiles_LP = 16#28#self.dataOut.nIncohInt_LP
5323 self.dataOut.nProfiles_LP = self.dataOut.data_acf.shape[1]#28#self.dataOut.nIncohInt_LP
5324 self.dataOut.NSCAN = 128
5325 self.dataOut.nIncohInt_LP = self.dataOut.nIncohInt*self.dataOut.NSCAN
5326 #print("sahpi",self.dataOut.nIncohInt_LP)
5327 #exit(1)
5328 self.dataOut.NLAG = 16
5329 self.dataOut.NRANGE = self.dataOut.data_acf.shape[-1]
5330
5331 #print(numpy.shape(self.dataOut.data_spc))
5332
5333 #exit(1)
5334 if mode==5:
5335 data = numpy.concatenate([getattr(data, attr_data) for data in data_inputs])
5336 setattr(self.dataOut, attr_data, data)
5337 data = numpy.concatenate([getattr(data, attr_data_2) for data in data_inputs])
5338 setattr(self.dataOut, attr_data_2, data)
5339 #data = numpy.concatenate([getattr(data, attr_data_3) for data in data_inputs])
5340 #setattr(self.dataOut, attr_data_3, data)
5341 #print(self.dataOut.moments.shape,self.dataOut.data_snr.shape,self.dataOut.heightList.shape) No newline at end of file
General Comments 0
You need to be logged in to leave comments. Login now