##// END OF EJS Templates
schain-jc
RSS
Noise.py...
Victor Sarmiento -
r83:e8a70c7ad5e4
parent child
Show More
Show file before | Show file after | Hide comments
@@ -0,0 +1,108
1 import numpy
2
3 from Model.JROHeader import *
4 from Model.Voltage import Voltage
5
6 def hildebrand_sekhon(Data, navg=1 ):
7 """
8 This method is for the objective determination of de noise level in Doppler spectra. This
9 implementation technique is based on the fact that the standard deviation of the spectral
10 densities is equal to the mean spectral density for white Gaussian noise
11
12 Inputs:
13 Data : heights
14 navg : numbers of averages
15
16 Return:
17 -1 : any error
18 anoise : noise's level
19 """
20 divisor = 8
21 ratio = 7 / divisor
22 data = Data.reshape(-1)
23 npts = data.size #numbers of points of the data
24
25 if npts < 32:
26 print "error in noise - requires at least 32 points"
27 return -1.0
28
29 # data sorted in ascending order
30 nmin = int(npts/divisor + ratio);
31 s = 0.0
32 s2 = 0.0
33 data2 = data[:npts]
34 data2.sort()
35
36 for i in range(nmin):
37 s += data2[i]
38 s2 += data2[i]**2;
39
40 icount = nmin
41 iflag = 0
42
43 for i in range(nmin, npts):
44 s += data2[i];
45 s2 += data2[i]**2
46 icount=icount+1;
47 p = s / float(icount);
48 p2 = p**2;
49 q = s2 / float(icount) - p2;
50 leftc = p2;
51 rightc = q * float(navg);
52
53 if leftc > rightc:
54 iflag = 1; #No weather signal
55 # Signal detect: R2 < 1 (R2 = leftc/rightc)
56 if(leftc < rightc):
57 if iflag:
58 break
59
60 anoise = 0.0;
61 for j in range(i):
62 anoise += data2[j];
63
64 anoise = anoise / float(i);
65
66 return anoise;
67
68
69 class Noise():
70 """
71 Clase que implementa los metodos necesarios para deternimar el nivel de ruido en un Spectro Doppler
72 """
73 m_DataObj = None
74
75
76 def __init__(self, m_Voltage=None):
77 """
78 Inicializador de la clase Noise para la la determinacion del nivel de ruido en un Spectro Doppler.
79
80 Affected:
81 self.m_DataObj
82
83 Return:
84 None
85 """
86 # if m_Voltage == None:
87 # m_Voltage = Voltage()
88 #
89 # if not(isinstance(m_Voltage, Voltage)):
90 # raise ValueError, "in Noise class, m_Voltage must be an Voltage class object"
91
92 self.m_DataObj = m_Voltage
93
94
95 def getNoiseLevelByHildebrandSekhon(self):
96 """
97 Determino el nivel de ruido usando el metodo Hildebrand-Sekhon
98
99 Return:
100 noise level
101 """
102 data = self.m_DataObj.data # heights x perfiles
103 #heights = numpy.transpose( Data, (2,0,1) ) # channel x profile x height
104 #data = Data[0,0,:]
105
106 noiselevel = hildebrand_sekhon(data)
107
108 print noiselevel
General Comments 0
You need to be logged in to leave comments. Login now