COMPASS  5.0.0
End-to-end AO simulation tool using GPU acceleration
correlation_bokeh Namespace Reference

Functions

def dphi_highpass (r, x0, tabx, taby)
 
def dphi_lowpass (r, x0, L0, tabx, taby)
 
def Ij0t83 (x, tabx, taby)
 
def unMoinsJ0 (x)
 
def tabulateIj0 (L0)
 
def asymp_macdo (x)
 
def macdo (x)
 
def rodconan (r, L0)
 
def variance (f, contributors, method="Default")
 Return the error variance of specified contributors params: f (h5py.File) : roket hdf5 file opened with h5py contributors (list of string) : list of the contributors method (optional, default="Default") : if "Independence", the function returns ths sum of the contributors variances. More...
 
def varianceMultiFiles (fs, frac_per_layer, contributors)
 Return the variance computed from the sum of contributors of roket files fs, ponderated by frac params: fs (list) : list of hdf5 files opened with h5py frac_per_layer (dict) : frac for each layer contributors (list of string) : list of the contributors return: v (np.array(dim=1)) : variance vector. More...
 
def cumulativeSR (v, Lambda_tar)
 Returns the cumulative Strehl ratio over the modes from the variance on each mode params: v (np.array(dim=1)) : variance vector return: s (np.array(dim=1)) : cumulative SR. More...
 
def update (attrs, old, new)
 

Variables

string datapath = "/home/fferreira/Data/correlation/"
 
 filenames = glob.glob(datapath + "roket_8m_1layer_dir*_cpu.h5")
 
list files = []
 
 ff = h5py.File(f, mode='r')
 
tuple nmodes = (files[0])["P"][:].shape[0]
 
list xpos = files[0].attrs["wfs.xpos"][0]
 
list ypos = files[0].attrs["wfs.ypos"][0]
 
list contributors = ["tomography", "bandwidth"]
 
list Lambda_tar = files[0].attrs["target.Lambda"][0]
 
list Lambda_wfs = files[0].attrs["wfs.Lambda"][0]
 
list L0 = files[0].attrs["L0"][0]
 
list dt = files[0].attrs["ittime"]
 
list H = files[0].attrs["atm.alt"][0]
 
int RASC = 180 / np.pi * 3600.
 
int Htheta
 
list r0 = files[0].attrs["r0"] * (Lambda_tar / Lambda_wfs)**(6. / 5.)
 
 nfiles = len(files)
 
 data = np.zeros((nmodes, 4, nfiles))
 
 theta = np.zeros(nfiles)
 
 speeds = np.zeros(nfiles)
 
 gain = np.zeros(nfiles)
 
 tabx
 
 taby
 
int ind = 0
 
 f
 
 method
 
tuple covar = (data[:, 0, :] - data[:, 3, :]) / 2.
 
 xaxis_select
 
 yaxis_select
 
 speed_select
 
 dir_select
 
 gain_select
 
 source = ColumnDataSource(data=dict(x=[], y=[], speed=[], theta=[], gain=[]))
 
 source_model = ColumnDataSource(data=dict(x=[], y=[], speed=[], theta=[], gain=[]))
 
 hover
 
string TOOLS = "resize,save,pan,box_zoom,tap, box_select, wheel_zoom, lasso_select,reset"
 
 p = figure(plot_height=600, plot_width=700, title="", tools=[hover, TOOLS])
 
 x
 
 y
 
 size
 
 color
 
dictionary xmap = {"Windspeed": speeds, "Winddir": theta, "Gain": gain}
 
dictionary ymap
 
list buttons = [xaxis_select, speed_select, dir_select, yaxis_select, gain_select]
 

Function Documentation

◆ asymp_macdo()

def correlation_bokeh.asymp_macdo (   x)

Definition at line 65 of file correlation_bokeh.py.

65 def asymp_macdo(x):
66  k2 = 1.00563491799858928388289314170833
67  k3 = 1.25331413731550012081
68  a1 = 0.22222222222222222222
69  a2 = -0.08641975308641974829
70  a3 = 0.08001828989483310284
71 
72  x_1 = 1. / x
73  res = k2 - k3 * np.exp(-x) * x**(1. / 3.) * (1.0 + x_1 * (a1 + x_1 *
74  (a2 + x_1 * a3)))
75  return res
76 
77 
Here is the caller graph for this function:

◆ cumulativeSR()

def correlation_bokeh.cumulativeSR (   v,
  Lambda_tar 
)

Returns the cumulative Strehl ratio over the modes from the variance on each mode params: v (np.array(dim=1)) : variance vector return: s (np.array(dim=1)) : cumulative SR.

Definition at line 190 of file correlation_bokeh.py.

190  """
191  s = np.cumsum(v)
192  s = np.exp(-s * (2 * np.pi / Lambda_tar)**2)
193 
194  return s
195 
196 

◆ dphi_highpass()

def correlation_bokeh.dphi_highpass (   r,
  x0,
  tabx,
  taby 
)

Definition at line 23 of file correlation_bokeh.py.

23 def dphi_highpass(r, x0, tabx, taby):
24  return (r**
25  (5. / 3.)) * (1.1183343328701949 - Ij0t83(r * (np.pi / x0), tabx, taby)) * (
26  2 * (2 * np.pi)**(8 / 3.) * 0.0228956)
27 
28 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ dphi_lowpass()

def correlation_bokeh.dphi_lowpass (   r,
  x0,
  L0,
  tabx,
  taby 
)

Definition at line 29 of file correlation_bokeh.py.

29 def dphi_lowpass(r, x0, L0, tabx, taby):
30  return rodconan(r, L0) - dphi_highpass(r, x0, tabx, taby)
31  #return (r**(5./3.)) * Ij0t83(r*(np.pi/x0), tabx, taby) * (2*(2*np.pi)**(8/3.)*0.0228956)
32 
33 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Ij0t83()

def correlation_bokeh.Ij0t83 (   x,
  tabx,
  taby 
)

Definition at line 34 of file correlation_bokeh.py.

34 def Ij0t83(x, tabx, taby):
35  if (x < np.exp(-3.0)):
36  return 0.75 * x**(1. / 3) * (1 - x**2 / 112.)
37  else:
38  return np.interp(x, tabx, taby)
39 
40 
Here is the caller graph for this function:

◆ macdo()

def correlation_bokeh.macdo (   x)

Definition at line 78 of file correlation_bokeh.py.

78 def macdo(x):
79  a = 5. / 6.
80  x2a = x**(2. * a)
81  x22 = x * x / 4.
82  s = 0.0
83 
84  Ga = [
85  0, 12.067619015983075, 5.17183672113560444, 0.795667187867016068,
86  0.0628158306210802181, 0.00301515986981185091, 9.72632216068338833e-05,
87  2.25320204494595251e-06, 3.93000356676612095e-08, 5.34694362825451923e-10,
88  5.83302941264329804e-12
89  ]
90 
91  Gma = [
92  -3.74878707653729304, -2.04479295083852408, -0.360845814853857083,
93  -0.0313778969438136685, -0.001622994669507603, -5.56455315259749673e-05,
94  -1.35720808599938951e-06, -2.47515152461894642e-08, -3.50257291219662472e-10,
95  -3.95770950530691961e-12, -3.65327031259100284e-14
96  ]
97 
98  x2n = 0.5
99 
100  s = Gma[0] * x2a
101  s *= x2n
102 
103  x2n *= x22
104 
105  for n in np.arange(10) + 1:
106 
107  s += (Gma[n] * x2a + Ga[n]) * x2n
108  x2n *= x22
109 
110  return s
111 
112 
Here is the caller graph for this function:

◆ rodconan()

def correlation_bokeh.rodconan (   r,
  L0 
)

Definition at line 113 of file correlation_bokeh.py.

113 def rodconan(r, L0):
114  res = 0
115  k1 = 0.1716613621245709486
116  dprf0 = (2 * np.pi / L0) * r
117  if (dprf0 > 4.71239):
118  res = asymp_macdo(dprf0)
119  else:
120  res = -macdo(dprf0)
121 
122  res *= k1 * L0**(5. / 3.)
123 
124  return res
125 
126 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ tabulateIj0()

def correlation_bokeh.tabulateIj0 (   L0)

Definition at line 49 of file correlation_bokeh.py.

49 def tabulateIj0(L0):
50  n = 10000
51  t = np.linspace(-4, 10, n)
52  dt = (t[-1] - t[0]) / (n - 1)
53  smallx = np.exp(-4.0)
54  A = 0.75 * smallx**(1. / 3) * (1 - smallx**2 / 112.)
55  X = np.exp(t)
56  #Y = np.exp(-t*(5./3.))*unMoinsJ0(X)
57  Y = (np.exp(2 * t) + (1. / L0)**2)**(-8. / 6.) * unMoinsJ0(X) * np.exp(t)
58  Y[1:] = np.cumsum(Y[:-1] + np.diff(Y) / 2.)
59  Y[0] = 0.
60  Y = Y * dt + A
61 
62  return X, Y
63 
64 
Here is the call graph for this function:

◆ unMoinsJ0()

def correlation_bokeh.unMoinsJ0 (   x)

Definition at line 41 of file correlation_bokeh.py.

41 def unMoinsJ0(x):
42  # if(x<0.1):
43  # x22 = (x/2.)**2
44  # return (1-x22/4.)*x22
45  # else:
46  return 1 - jv(0, x)
47 
48 
Here is the caller graph for this function:

◆ update()

def correlation_bokeh.update (   attrs,
  old,
  new 
)

Definition at line 197 of file correlation_bokeh.py.

197 def update(attrs, old, new):
198  speed = speed_select.value
199  direction = dir_select.value
200  g = gain_select.value
201  xname = xaxis_select.value
202  yname = yaxis_select.value
203 
204  ydata = ymap[yname]
205  x = xmap[xname]
206 
207  ind = np.ones(ydata.shape[0])
208  if (direction != "All"):
209  ind *= (xmap["Winddir"] == float(direction))
210  if (speed != "All"):
211  ind *= (xmap["Windspeed"] == float(speed))
212  if (g != "All"):
213  ind *= (xmap["Gain"] == float(g))
214 
215  ind = np.where(ind)
216  if (yname == b"Var(t)"):
217  Hthetak = Htheta / xmap["Gain"]
218  y_model = np.ones(ind[0].shape[0])
219  #y_model = y_model * 6.88 * (Htheta/r0)**(5./3.) * 0.5
220  for k in range(ind[0].shape[0]):
221  y_model[k] = dphi_lowpass(Htheta, 0.2, L0, tabx, taby) * (1 / r0)**(
222  5. / 3.) * 0.5 #* xmap["Gain"][ind][k]**2
223  if (yname == b"Var(bp)"):
224  vdt = xmap["Windspeed"] * dt / xmap["Gain"]
225  y_model = np.zeros(vdt[ind].shape[0])
226  for k in range(vdt[ind].shape[0]):
227  y_model[k] = dphi_lowpass(vdt[ind][k], 0.2, L0, tabx,
228  taby) * (1. / r0)**(5. / 3.) * 0.5
229  if (yname == b"Covar"):
230  vdt = xmap["Windspeed"] * dt / xmap["Gain"]
231  Hthetak = Htheta / xmap["Gain"]
232  gamma = np.arctan2(ypos, xpos) - xmap["Winddir"] * np.pi / 180.
233  rho = np.sqrt(Htheta**2 + (vdt)**2 - 2 * Htheta * vdt * np.cos(gamma))
234  Drho = np.zeros(rho[ind].shape[0])
235  Dt = Drho.copy()
236  for k in range(rho[ind].shape[0]):
237  Drho[k] = dphi_lowpass(rho[ind][k], 0.2, L0, tabx,
238  taby) * (1 / r0)**(5. / 3.)
239  #Drho = 6.88 * (rho[ind]/r0)**(5./3.)
240  for k in range(Dt.shape[0]):
241  Dt[k] = dphi_lowpass(Htheta, 0.2, L0, tabx, taby) * (1 / r0)**(
242  5. / 3.) # * xmap["Gain"][ind][k]**2
243  #Dt = 6.88 * (Htheta/r0)**(5./3.)
244  Dbp = np.zeros(vdt[ind].shape[0])
245  for k in range(vdt[ind].shape[0]):
246  Dbp[k] = dphi_lowpass(vdt[ind][k], 0.2, L0, tabx, taby) * (1 / r0)**(5. / 3.)
247  #Dbp = 6.88 * (vdt[ind]/r0) ** (5./3.)
248  y_model = 0.5 * (Dt + Dbp - Drho)
249 
250  source.data = dict(x=x[ind], y=ydata[ind], speed=xmap["Windspeed"][ind],
251  theta=xmap["Winddir"][ind], gain=xmap["Gain"][ind])
252  source_model.data = dict(x=x[ind], y=y_model, speed=xmap["Windspeed"][ind],
253  theta=xmap["Winddir"][ind], gain=xmap["Gain"][ind])
254 
255 
Here is the call graph for this function:

◆ variance()

def correlation_bokeh.variance (   f,
  contributors,
  method = "Default" 
)

Return the error variance of specified contributors params: f (h5py.File) : roket hdf5 file opened with h5py contributors (list of string) : list of the contributors method (optional, default="Default") : if "Independence", the function returns ths sum of the contributors variances.

If "Default", it returns the variance of the contributors sum

Definition at line 135 of file correlation_bokeh.py.

135  """
136  P = f["P"][:]
137  nmodes = P.shape[0]
138  swap = np.arange(nmodes) - 2
139  swap[0:2] = [nmodes - 2, nmodes - 1]
140  if (method == b"Default"):
141  err = f[contributors[0]][:] * 0.
142  for c in contributors:
143  err += f[c][:]
144  return np.var(P.dot(err), axis=1)[swap], np.var(
145  P.dot(f["tomography"][:]),
146  axis=1)[swap], np.var(P.dot(f["bandwidth"][:]), axis=1)[swap]
147 
148  elif (method == b"Independence"):
149  nmodes = P.shape[0]
150  v = np.zeros(nmodes)
151  for c in contributors:
152  v += np.var(P.dot(f[c][:]), axis=1)
153  return v[swap]
154 
155  else:
156  raise TypeError("Wrong method input")
157 
158 

◆ varianceMultiFiles()

def correlation_bokeh.varianceMultiFiles (   fs,
  frac_per_layer,
  contributors 
)

Return the variance computed from the sum of contributors of roket files fs, ponderated by frac params: fs (list) : list of hdf5 files opened with h5py frac_per_layer (dict) : frac for each layer contributors (list of string) : list of the contributors return: v (np.array(dim=1)) : variance vector.

Definition at line 168 of file correlation_bokeh.py.

168  """
169  f = fs[0]
170  P = f["P"][:]
171  nmodes = P.shape[0]
172  swap = np.arange(nmodes) - 2
173  swap[0:2] = [nmodes - 2, nmodes - 1]
174  err = f[contributors[0]][:] * 0.
175  for f in fs:
176  frac = frac_per_layer[f.attrs["atm.alt"][0]]
177  for c in contributors:
178  err += np.sqrt(frac) * f[c][:]
179 
180  return np.var(P.dot(err), axis=1)[swap]
181 
182 

Variable Documentation

◆ buttons

list correlation_bokeh.buttons = [xaxis_select, speed_select, dir_select, yaxis_select, gain_select]

Definition at line 333 of file correlation_bokeh.py.

◆ color

correlation_bokeh.color

Definition at line 321 of file correlation_bokeh.py.

◆ contributors

correlation_bokeh.contributors = ["tomography", "bandwidth"]

Definition at line 268 of file correlation_bokeh.py.

◆ covar

tuple correlation_bokeh.covar = (data[:, 0, :] - data[:, 3, :]) / 2.

Definition at line 300 of file correlation_bokeh.py.

◆ data

correlation_bokeh.data = np.zeros((nmodes, 4, nfiles))

Definition at line 280 of file correlation_bokeh.py.

◆ datapath

string correlation_bokeh.datapath = "/home/fferreira/Data/correlation/"

Definition at line 256 of file correlation_bokeh.py.

◆ dir_select

correlation_bokeh.dir_select
Initial value:
1 = Select(title="Winddirs", value="All",
2  options=["All"] + [str(s) for s in np.unique(theta)])

Definition at line 310 of file correlation_bokeh.py.

◆ dt

list correlation_bokeh.dt = files[0].attrs["ittime"]

Definition at line 272 of file correlation_bokeh.py.

◆ f

correlation_bokeh.f

Definition at line 294 of file correlation_bokeh.py.

◆ ff

correlation_bokeh.ff = h5py.File(f, mode='r')

Definition at line 261 of file correlation_bokeh.py.

◆ filenames

correlation_bokeh.filenames = glob.glob(datapath + "roket_8m_1layer_dir*_cpu.h5")

Definition at line 257 of file correlation_bokeh.py.

◆ files

list correlation_bokeh.files = []

Definition at line 259 of file correlation_bokeh.py.

◆ gain

correlation_bokeh.gain = np.zeros(nfiles)

Definition at line 283 of file correlation_bokeh.py.

◆ gain_select

correlation_bokeh.gain_select
Initial value:
1 = Select(title="Gain", value="All",
2  options=["All"] + [str(s)[:3] for s in np.unique(gain)])

Definition at line 312 of file correlation_bokeh.py.

◆ H

list correlation_bokeh.H = files[0].attrs["atm.alt"][0]

Definition at line 273 of file correlation_bokeh.py.

◆ hover

correlation_bokeh.hover
Initial value:
1 = HoverTool(tooltips=[("Speed", "@speed"), ("Winddir", "@theta"), ("Gain",
2  "@gain")])

Definition at line 316 of file correlation_bokeh.py.

◆ Htheta

int correlation_bokeh.Htheta
Initial value:
1 = np.linalg.norm(
2  [xpos, ypos]
3 ) / RASC * H

Definition at line 275 of file correlation_bokeh.py.

◆ ind

int correlation_bokeh.ind = 0

Definition at line 291 of file correlation_bokeh.py.

◆ L0

list correlation_bokeh.L0 = files[0].attrs["L0"][0]

Definition at line 271 of file correlation_bokeh.py.

◆ Lambda_tar

list correlation_bokeh.Lambda_tar = files[0].attrs["target.Lambda"][0]

Definition at line 269 of file correlation_bokeh.py.

◆ Lambda_wfs

list correlation_bokeh.Lambda_wfs = files[0].attrs["wfs.Lambda"][0]

Definition at line 270 of file correlation_bokeh.py.

◆ method

correlation_bokeh.method

Definition at line 294 of file correlation_bokeh.py.

◆ nfiles

correlation_bokeh.nfiles = len(files)

Definition at line 279 of file correlation_bokeh.py.

◆ nmodes

tuple correlation_bokeh.nmodes = (files[0])["P"][:].shape[0]

Definition at line 265 of file correlation_bokeh.py.

◆ p

correlation_bokeh.p = figure(plot_height=600, plot_width=700, title="", tools=[hover, TOOLS])

Definition at line 320 of file correlation_bokeh.py.

◆ r0

list correlation_bokeh.r0 = files[0].attrs["r0"] * (Lambda_tar / Lambda_wfs)**(6. / 5.)

Definition at line 278 of file correlation_bokeh.py.

◆ RASC

int correlation_bokeh.RASC = 180 / np.pi * 3600.

Definition at line 274 of file correlation_bokeh.py.

◆ size

correlation_bokeh.size

Definition at line 321 of file correlation_bokeh.py.

◆ source

correlation_bokeh.source = ColumnDataSource(data=dict(x=[], y=[], speed=[], theta=[], gain=[]))

Definition at line 314 of file correlation_bokeh.py.

◆ source_model

correlation_bokeh.source_model = ColumnDataSource(data=dict(x=[], y=[], speed=[], theta=[], gain=[]))

Definition at line 315 of file correlation_bokeh.py.

◆ speed_select

correlation_bokeh.speed_select
Initial value:
1 = Select(title="Windspeeds", value="All",
2  options=["All"] + [str(s) for s in np.unique(speeds)])

Definition at line 308 of file correlation_bokeh.py.

◆ speeds

correlation_bokeh.speeds = np.zeros(nfiles)

Definition at line 282 of file correlation_bokeh.py.

◆ tabx

correlation_bokeh.tabx

Definition at line 285 of file correlation_bokeh.py.

◆ taby

correlation_bokeh.taby

Definition at line 285 of file correlation_bokeh.py.

◆ theta

correlation_bokeh.theta = np.zeros(nfiles)

Definition at line 281 of file correlation_bokeh.py.

◆ TOOLS

string correlation_bokeh.TOOLS = "resize,save,pan,box_zoom,tap, box_select, wheel_zoom, lasso_select,reset"

Definition at line 318 of file correlation_bokeh.py.

◆ x

correlation_bokeh.x

Definition at line 321 of file correlation_bokeh.py.

◆ xaxis_select

correlation_bokeh.xaxis_select
Initial value:
1 = Select(title="X-axis", value="Windspeed",
2  options=["Windspeed", "Winddir", "Gain"])

Definition at line 302 of file correlation_bokeh.py.

◆ xmap

dictionary correlation_bokeh.xmap = {"Windspeed": speeds, "Winddir": theta, "Gain": gain}

Definition at line 324 of file correlation_bokeh.py.

◆ xpos

list correlation_bokeh.xpos = files[0].attrs["wfs.xpos"][0]

Definition at line 266 of file correlation_bokeh.py.

◆ y

correlation_bokeh.y

Definition at line 321 of file correlation_bokeh.py.

◆ yaxis_select

correlation_bokeh.yaxis_select
Initial value:
1 = Select(
2  title="Y-axis", value="Covar",
3  options=["Covar", "Var(t+bp)", "Var(t)", "Var(bp)", "Var(t)+Var(bp)"])

Definition at line 304 of file correlation_bokeh.py.

◆ ymap

dictionary correlation_bokeh.ymap
Initial value:
1 = {
2  "Covar": np.sum(covar, axis=0),
3  "Var(t+bp)": np.sum(data[:, 0, :], axis=0),
4  "Var(t)": np.sum(data[:, 1, :], axis=0),
5  "Var(bp)": np.sum(data[:, 2, :], axis=0),
6  "Var(t)+Var(bp)": np.sum(data[:, 3, :], axis=0)
7 }

Definition at line 325 of file correlation_bokeh.py.

◆ ypos

list correlation_bokeh.ypos = files[0].attrs["wfs.ypos"][0]

Definition at line 267 of file correlation_bokeh.py.