| @@ -665,11 +665,13 class Plot(Operation): | |||
|
|
665 | 665 | tm = getattr(dataOut, self.attr_time) |
|
|
666 | 666 | |
|
|
667 | 667 | if self.data and 'time' in self.xaxis and (tm - self.tmin) >= self.xrange*60*60: |
|
|
668 | self.clear_figures() | |
|
|
668 | 669 | self.save_time = tm |
|
|
670 | self.__plot() | |
|
|
669 | 671 | self.tmin += self.xrange*60*60 |
|
|
670 | 672 | self.data.setup() |
|
|
671 | self.clear_figures() | |
|
|
672 |
|
|
|
|
673 | #self.clear_figures() | |
|
|
674 | ||
|
|
673 | 675 | |
|
|
674 | 676 | self.__update(dataOut, tm) |
|
|
675 | 677 | |
| @@ -320,36 +320,47 class selectHeights(Operation): | |||
|
|
320 | 320 | |
|
|
321 | 321 | |
|
|
322 | 322 | class filterByHeights(Operation): |
|
|
323 | ||
|
|
323 | ifConfig=False | |
|
|
324 | deltaHeight = None | |
|
|
325 | newdelta=None | |
|
|
326 | newheights=None | |
|
|
327 | r=None | |
|
|
328 | h0=None | |
|
|
329 | nHeights=None | |
|
|
324 | 330 | def run(self, dataOut, window): |
|
|
325 | 331 | |
|
|
326 | deltaHeight = dataOut.heightList[1] - dataOut.heightList[0] | |
|
|
327 | 332 | |
|
|
333 | # print("1",dataOut.data.shape) | |
|
|
334 | # print(dataOut.nHeights) | |
|
|
328 | 335 | if window == None: |
|
|
329 | window = (dataOut.radarControllerHeaderObj.txA/dataOut.radarControllerHeaderObj.nBaud) / deltaHeight | |
|
|
330 | ||
|
|
331 | newdelta = deltaHeight * window | |
|
|
332 | r = dataOut.nHeights % window | |
|
|
333 | newheights = (dataOut.nHeights-r)/window | |
|
|
336 | window = (dataOut.radarControllerHeaderObj.txA/dataOut.radarControllerHeaderObj.nBaud) / self.deltaHeight | |
|
|
334 | 337 | |
|
|
335 | if newheights <= 1: | |
|
|
338 | if not self.ifConfig: #and dataOut.useInputBuffer: | |
|
|
339 | self.deltaHeight = dataOut.heightList[1] - dataOut.heightList[0] | |
|
|
340 | self.ifConfig = True | |
|
|
341 | self.newdelta = self.deltaHeight * window | |
|
|
342 | self.r = dataOut.nHeights % window | |
|
|
343 | self.newheights = (dataOut.nHeights-self.r)/window | |
|
|
344 | self.h0 = dataOut.heightList[0] | |
|
|
345 | self.nHeights = dataOut.nHeights | |
|
|
346 | if self.newheights <= 1: | |
|
|
336 | 347 | raise ValueError("filterByHeights: Too few heights. Current number of heights is %d and window is %d" %(dataOut.nHeights, window)) |
|
|
337 | 348 | |
|
|
338 | 349 | if dataOut.flagDataAsBlock: |
|
|
339 | 350 | """ |
|
|
340 | 351 | Si la data es obtenida por bloques, dimension = [nChannels, nProfiles, nHeis] |
|
|
341 | 352 | """ |
|
|
342 |
buffer = dataOut.data[:, :, 0:int( |
|
|
|
343 |
buffer = buffer.reshape(dataOut.nChannels, dataOut.nProfiles, int( |
|
|
|
353 | buffer = dataOut.data[:, :, 0:int(self.nHeights-self.r)] | |
|
|
354 | buffer = buffer.reshape(dataOut.nChannels, dataOut.nProfiles, int(self.nHeights/window), window) | |
|
|
344 | 355 | buffer = numpy.sum(buffer,3) |
|
|
345 | 356 | |
|
|
346 | 357 | else: |
|
|
347 |
buffer = dataOut.data[:,0:int( |
|
|
|
348 |
buffer = buffer.reshape(dataOut.nChannels,int( |
|
|
|
358 | buffer = dataOut.data[:,0:int(self.nHeights-self.r)] | |
|
|
359 | buffer = buffer.reshape(dataOut.nChannels,int(self.nHeights/window),int(window)) | |
|
|
349 | 360 | buffer = numpy.sum(buffer,2) |
|
|
350 | 361 | |
|
|
351 | 362 | dataOut.data = buffer |
|
|
352 |
dataOut.heightList = |
|
|
|
363 | dataOut.heightList = self.h0 + numpy.arange( self.newheights )*self.newdelta | |
|
|
353 | 364 | dataOut.windowOfFilter = window |
|
|
354 | 365 | |
|
|
355 | 366 | #update Processing Header: |
| @@ -359,6 +370,7 class filterByHeights(Operation): | |||
|
|
359 | 370 | return dataOut |
|
|
360 | 371 | |
|
|
361 | 372 | |
|
|
373 | ||
|
|
362 | 374 | class setH0(Operation): |
|
|
363 | 375 | |
|
|
364 | 376 | def run(self, dataOut, h0, deltaHeight = None): |
| @@ -3276,7 +3288,7 class RemoveProfileSats2(Operation): | |||
|
|
3276 | 3288 | |
|
|
3277 | 3289 | |
|
|
3278 | 3290 | __slots__ = ('n','navg','profileMargin','thHistOutlier','minHei_idx','maxHei_idx','nHeights', |
|
|
3279 | 'first_utcBlock','__profIndex','init_prof','end_prof','lenProfileOut','nChannels', | |
|
|
3291 | 'first_utcBlock','__profIndex','init_prof','end_prof','lenProfileOut','nChannels','cohFactor', | |
|
|
3280 | 3292 | '__count_exec','__initime','__dataReady','__ipp', 'minRef', 'maxRef', 'debug','prev_pnoise') |
|
|
3281 | 3293 | def __init__(self, **kwargs): |
|
|
3282 | 3294 | |
| @@ -3284,7 +3296,7 class RemoveProfileSats2(Operation): | |||
|
|
3284 | 3296 | self.isConfig = False |
|
|
3285 | 3297 | self.currentTime = None |
|
|
3286 | 3298 | |
|
|
3287 |
def setup(self,dataOut, n=None , navg=0. |
|
|
|
3299 | def setup(self,dataOut, n=None , navg=0.9, profileMargin=50,thHistOutlier=15,minHei=None, maxHei=None, nBins=10, | |
|
|
3288 | 3300 | minRef=None, maxRef=None, debug=False, remYagi=False, nProfYagi = 0, offYagi=0, minHJULIA=None, maxHJULIA=None, |
|
|
3289 | 3301 | idate=None,startH=None,endH=None ): |
|
|
3290 | 3302 | |
| @@ -3323,7 +3335,7 class RemoveProfileSats2(Operation): | |||
|
|
3323 | 3335 | self.test_counter = 0 |
|
|
3324 | 3336 | self.debug = debug |
|
|
3325 | 3337 | self.remYagi = remYagi |
|
|
3326 | ||
|
|
3338 | self.cohFactor = dataOut.nCohInt | |
|
|
3327 | 3339 | if self.remYagi : |
|
|
3328 | 3340 | if minHJULIA==None or maxHJULIA==None: |
|
|
3329 | 3341 | raise ValueError("Parameters minHYagi and minHYagi are necessary!") |
| @@ -3338,7 +3350,7 class RemoveProfileSats2(Operation): | |||
|
|
3338 | 3350 | self.startTime = datetime.datetime.combine(idate,startH) |
|
|
3339 | 3351 | self.endTime = datetime.datetime.combine(idate,endH) |
|
|
3340 | 3352 | |
|
|
3341 |
log.warning("Be careful with the selection of parameters for sat removal! I |
|
|
|
3353 | log.warning("Be careful with the selection of parameters for sats removal! It is avisable to \ | |
|
|
3342 | 3354 | activate the debug parameter in this operation for calibration", self.name) |
|
|
3343 | 3355 | |
|
|
3344 | 3356 | |
| @@ -3370,14 +3382,19 activate the debug parameter in this operation for calibration", self.name) | |||
|
|
3370 | 3382 | #print(self.min_ref,self.max_ref) |
|
|
3371 | 3383 | |
|
|
3372 | 3384 | import scipy.signal |
|
|
3373 |
b, a = scipy.signal.butter(3, 0. |
|
|
|
3374 |
noise_ref = (data[c,:,self.min_ref:self.max_ref] * numpy.conjugate(data[c,:,self.min_ref:self.max_ref])) |
|
|
|
3385 | b, a = scipy.signal.butter(3, 0.5) | |
|
|
3386 | #noise_ref = (data[c,:,self.min_ref:self.max_ref] * numpy.conjugate(data[c,:,self.min_ref:self.max_ref])) | |
|
|
3387 | noise_ref = numpy.abs(data[c,:,self.min_ref:self.max_ref]) | |
|
|
3388 | lnoise = len(noise_ref[0,:]) | |
|
|
3389 | #print(noise_ref.shape) | |
|
|
3375 | 3390 | noise_ref = noise_ref.mean(axis=1) |
|
|
3391 | #fnoise = noise_ref | |
|
|
3376 | 3392 | fnoise = scipy.signal.filtfilt(b, a, noise_ref) |
|
|
3377 | 3393 | #noise_refdB = 10* numpy.log10(noise_ref) |
|
|
3378 | 3394 | #print("Noise ",numpy.percentile(noise_ref,95)) |
|
|
3379 |
p |
|
|
|
3395 | p95 = numpy.percentile(fnoise,90) | |
|
|
3380 | 3396 | mean_noise = fnoise.mean() |
|
|
3397 | ||
|
|
3381 | 3398 | if self.prev_pnoise != None: |
|
|
3382 | 3399 | if mean_noise < (1.5 * self.prev_pnoise) : |
|
|
3383 | 3400 | self.prev_pnoise = mean_noise |
| @@ -3390,20 +3407,24 activate the debug parameter in this operation for calibration", self.name) | |||
|
|
3390 | 3407 | |
|
|
3391 | 3408 | |
|
|
3392 | 3409 | |
|
|
3393 |
power = (data[c,:,self.minHei_idx:self.maxHei_idx] * numpy.conjugate(data[c,:,self.minHei_idx:self.maxHei_idx])) |
|
|
|
3394 | heis = len(power[0,:]) | |
|
|
3410 | #power = (data[c,:,self.minHei_idx:self.maxHei_idx] * numpy.conjugate(data[c,:,self.minHei_idx:self.maxHei_idx])) | |
|
|
3411 | power = numpy.abs(data[c,:,self.minHei_idx:self.maxHei_idx]) | |
|
|
3412 | npower = len(power[0,:]) | |
|
|
3413 | #print(power.shape) | |
|
|
3395 | 3414 | power = power.mean(axis=1) |
|
|
3396 | 3415 | |
|
|
3397 | 3416 | fpower = scipy.signal.filtfilt(b, a, power) |
|
|
3398 | 3417 | #print(power.shape) |
|
|
3399 | 3418 | #powerdB = 10* numpy.log10(power) |
|
|
3400 | 3419 | |
|
|
3401 |
th = p |
|
|
|
3420 | th = p95 #* 1.1 | |
|
|
3421 | #th = mean_noise | |
|
|
3402 | 3422 | index = numpy.where(fpower > th ) |
|
|
3403 |
#print("Noise ",mean_noise, p |
|
|
|
3423 | #print("Noise ",mean_noise, p95) | |
|
|
3404 | 3424 | #print(index) |
|
|
3405 | 3425 | |
|
|
3406 | if index[0].size < int(self.navg*profiles): | |
|
|
3426 | ||
|
|
3427 | if index[0].size <= int(self.navg*profiles): | |
|
|
3407 | 3428 | indexes = numpy.append(indexes, index[0]) |
|
|
3408 | 3429 | |
|
|
3409 | 3430 | #print("sdas ", noise_ref.mean()) |
| @@ -3423,12 +3444,15 activate the debug parameter in this operation for calibration", self.name) | |||
|
|
3423 | 3444 | # plt.grid() |
|
|
3424 | 3445 | # plt.show() |
|
|
3425 | 3446 | |
|
|
3426 | # fig, ax = plt.subplots() | |
|
|
3427 |
|
|
|
|
3428 | # ax.axhline(th, color='r') | |
|
|
3429 | # ax.axhline(std, color='b') | |
|
|
3430 | # plt.grid() | |
|
|
3431 | # plt.show() | |
|
|
3447 | if self.debug: | |
|
|
3448 | fig, ax = plt.subplots() | |
|
|
3449 | ax.plot(fpower, label="power") | |
|
|
3450 | #ax.plot(fnoise, label="noise ref") | |
|
|
3451 | ax.axhline(th, color='g', label="th") | |
|
|
3452 | #ax.axhline(std, color='b', label="mean") | |
|
|
3453 | ax.legend() | |
|
|
3454 | plt.grid() | |
|
|
3455 | plt.show() | |
|
|
3432 | 3456 | |
|
|
3433 | 3457 | #print(indexes) |
|
|
3434 | 3458 | |
| @@ -3464,16 +3488,18 activate the debug parameter in this operation for calibration", self.name) | |||
|
|
3464 | 3488 | indexesYagi = indexesYagi[ (indexesYagi >= 0) & (indexesYagi<profiles)] |
|
|
3465 | 3489 | indexesYagi = numpy.unique(indexesYagi) |
|
|
3466 | 3490 | |
|
|
3491 | #print("indexes: " ,indexes) | |
|
|
3467 | 3492 | outs_lines = numpy.unique(indexes) |
|
|
3468 | ||
|
|
3493 | #print(outs_lines) | |
|
|
3469 | 3494 | |
|
|
3470 | 3495 | #Agrupando el histograma de outliers, |
|
|
3471 | 3496 | my_bins = numpy.linspace(0,int(profiles), div, endpoint=True) |
|
|
3472 | 3497 | |
|
|
3473 | 3498 | |
|
|
3474 | 3499 | hist, bins = numpy.histogram(outs_lines,bins=my_bins) |
|
|
3475 | hist_outliers_indexes = numpy.where(hist >= self.thHistOutlier) #es outlier | |
|
|
3476 | hist_outliers_indexes = hist_outliers_indexes[0] | |
|
|
3500 | #print("hist: ",hist) | |
|
|
3501 | hist_outliers_indexes = numpy.where(hist >= self.thHistOutlier)[0] #es outlier | |
|
|
3502 | #print(hist_outliers_indexes) | |
|
|
3477 | 3503 | if len(hist_outliers_indexes>0): |
|
|
3478 | 3504 | hist_outliers_indexes = numpy.append(hist_outliers_indexes,hist_outliers_indexes[-1]+1) |
|
|
3479 | 3505 | |
| @@ -3516,7 +3542,7 activate the debug parameter in this operation for calibration", self.name) | |||
|
|
3516 | 3542 | if self.remYagi : |
|
|
3517 | 3543 | ax[0].vlines(indexesYagi,750,850, linestyles='dashed', label = 'yagi', color='m') |
|
|
3518 | 3544 | else: |
|
|
3519 | c = ax[0].pcolormesh(x, y, dat.T, cmap ='jet', vmin = 60, vmax = 70) | |
|
|
3545 | c = ax[0].pcolormesh(x, y, dat.T, cmap ='jet', vmin = 60, vmax = (70+2*self.cohFactor)) | |
|
|
3520 | 3546 | ax[0].vlines(outs_lines,650,700, linestyles='dashed', label = 'outs', color='w') |
|
|
3521 | 3547 | #fig.colorbar(c) |
|
|
3522 | 3548 | ax[0].vlines(outlier_loc_index,700,750, linestyles='dashed', label = 'outs', color='r') |
General Comments 0
You need to be logged in to leave comments.
Login now