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

Functions

def compute_psf (filename)
 
def compute_psf_independence (filename)
 
def compute_and_compare_PSFs (filename, plot=False)
 
def filter_piston_TT (filename, C)
 
def filter_TT (filename, C)
 
def compute_covariance_model (filename)
 
def add_TT_model (filename, Ccov)
 
def load_datas (files)
 
def ensquare_PSF (filename, psf, N, display=False)
 
def cutsPSF (filename, psf, psfs)
 
def Hcor (f, Fe, g, dt)
 
def Hretard (f, Fe, g, dt)
 
def compareTransferFunctions (filename)
 

Variables

dictionary font = {'family': 'normal', 'weight': 'bold', 'size': 22}
 
string datapath = '/home/fferreira/Data/correlation/'
 
 filenames = glob.glob(datapath + 'roket_8m_1layer_dir*_cpu.h5')
 
list files = []
 
 tabx
 
 taby
 
 nfiles = len(filenames)
 
 theta = np.zeros(nfiles)
 
 speeds = np.zeros(nfiles)
 
 gain = np.zeros(nfiles)
 
 SRcompass = np.zeros(nfiles)
 
 SRroket = np.zeros(nfiles)
 
 SRi = np.zeros(nfiles)
 
 fROKET = h5py.File('ROKETStudy.h5', 'r')
 
 psfr = fROKET["psf"][:]
 
 psfi = fROKET["psfi"][:]
 
 nrjcompass = np.zeros(nfiles)
 
 nrjroket = np.zeros(nfiles)
 
 nrji = np.zeros(nfiles)
 
int ind = 0
 
 f = h5py.File('corStudy_Nact.h5', 'r')
 
 psf = f["psf"][:]
 
 psfs = f["psfs"][:]
 
 nrj = f["nrj5"][:]
 
 nrjs = f["nrj5s"][:]
 
 SR = np.max(psf, axis=(0, 1))
 
 SRs = np.max(psfs, axis=(0, 1))
 
list colors = ["blue", "red", "green"]
 
int k = 0
 
list c = colors[i]
 
 v = np.unique(speeds)[i]
 
 color
 
 s
 
 loc
 

Function Documentation

◆ add_TT_model()

def correlation_study.add_TT_model (   filename,
  Ccov 
)

Definition at line 246 of file correlation_study.py.

246 def add_TT_model(filename, Ccov):
247  C = np.zeros((Ccov.shape[0] + 2, Ccov.shape[0] + 2))
248  IF, T = rexp.get_IF(filename)
249  IF = IF.T
250  T = T.T
251  N = IF.shape[0]
252  deltaTT = T.T.dot(T) / N
253  deltaF = IF.T.dot(T) / N
254  pzt2tt = np.linalg.inv(deltaTT).dot(deltaF.T)
255 
256  CTT = Ccov - filter_TT(filename, Ccov)
257  C[-2:, -2:] = pzt2tt.dot(CTT).dot(pzt2tt.T)
258  C[:-2, :-2] = Ccov
259 
260  return C
261 
262 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compareTransferFunctions()

def correlation_study.compareTransferFunctions (   filename)

Definition at line 373 of file correlation_study.py.

373 def compareTransferFunctions(filename):
374  rfile = h5py.File(filename, 'r')
375  v = rfile.attrs["windspeed"][0]
376  dt = rfile.attrs["ittime"]
377  Fe = 1 / dt
378  g = rfile.attrs["gain"]
379  Lambda_tar = rfile.attrs["target.Lambda"][0]
380  Lambda_wfs = rfile.attrs["wfs.Lambda"][0]
381  r0 = rfile.attrs["r0"] * (Lambda_tar / Lambda_wfs)**(6. / 5.)
382  d = rfile.attrs["tel_diam"] / rfile.attrs["nxsub"]
383  fc = 0.314 * v / d
384  f = np.linspace(0.1, fc * 1.5, 1000)
385  H = Hcor(f, Fe, g, dt)
386  Hr = Hretard(f, Fe, g, dt)
387  plt.figure()
388  plt.plot(f, H)
389  plt.plot(f, Hr, color="red")
390  plt.plot([fc, fc], [H.min(), Hr.max()])
391  rfile.close()
392 
393 
Here is the call graph for this function:

◆ compute_and_compare_PSFs()

def correlation_study.compute_and_compare_PSFs (   filename,
  plot = False 
)

Definition at line 37 of file correlation_study.py.

37 def compute_and_compare_PSFs(filename, plot=False):
38  f = h5py.File(filename, 'r')
39  psf_compass = np.fft.fftshift(f["psf"][:])
40  #psf = compute_psf(filename)
41  #psfi = compute_psf_independence(filename)
42  psfc = 0
43  psfs = 0
44  tic = time.time()
45  Caniso, Cbp, Ccov = compute_covariance_model(filename)
46 
47  Nact = f["Nact"][:]
48  N1 = np.linalg.inv(Nact)
49  Caniso = N1.dot(Caniso).dot(N1)
50  Cbp = N1.dot(Cbp).dot(N1)
51  Ccov = N1.dot(Ccov).dot(N1)
52 
53  Ccov_filtered = filter_piston_TT(filename, Ccov)
54  Ctt = add_TT_model(filename, Ccov)
55  Ctt[:-2, :-2] = Ccov_filtered
56  Ctt = Ctt + Ctt.T
57  Caniso_filtered = filter_piston_TT(filename, Caniso)
58  tmp = add_TT_model(filename, Caniso)
59  tmp[:-2, :-2] = Caniso_filtered
60  Ctt += tmp
61  Cbp_filtered = filter_piston_TT(filename, Cbp)
62  tmp = add_TT_model(filename, Cbp)
63  tmp[:-2, :-2] = Cbp_filtered
64  Ctt += tmp
65  #contributors = ["noise","aliasing","non linearity","filtered modes"]
66  #cov_err = rexp.get_coverr_independence_contributors(filename,contributors)
67  P = f["P"][:]
68  cov_err = P.dot(Ctt).dot(P.T)
69  otftels, otf2s, psfs, gpu = gamora.psf_rec_Vii(filename, fitting=False,
70  cov=cov_err.astype(np.float32))
71  tac = time.time()
72  print("PSF estimated in ", tac - tic, " seconds")
73  t = f["tomography"][:]
74  b = f["bandwidth"][:]
75  tb = t + b
76  tb = tb.dot(tb.T) / float(tb.shape[1])
77  cov_err = P.dot(tb).dot(P.T)
78  otftel, otf2, psf, gpu = gamora.psf_rec_Vii(filename, fitting=False,
79  cov=cov_err.astype(np.float32))
80  if (plot):
81  Lambda_tar = f.attrs["target.Lambda"][0]
82  RASC = 180 / np.pi * 3600.
83  pixsize = Lambda_tar * 1e-6 / (
84  psf.shape[0] * f.attrs["tel_diam"] / f.attrs["pupdiam"]) * RASC
85  x = (np.arange(psf.shape[0]) - psf.shape[0] / 2) * pixsize / (
86  Lambda_tar * 1e-6 / f.attrs["tel_diam"] * RASC)
87  plt.figure()
88  plt.subplot(2, 1, 1)
89  plt.semilogy(x, psf[psf.shape[0] / 2, :], color="blue")
90  plt.semilogy(x, psfs[psf.shape[0] / 2, :], color="red")
91  plt.xlabel("X-axis angular distance [units of lambda/D]")
92  plt.ylabel("Normalized intensity")
93  plt.legend(["PSF exp", "PSF model"])
94 
95  #plt.semilogy(x,psf_compass[psf.shape[0]/2,:],color="red")
96  plt.legend(["PSF exp", "PSF model"])
97  plt.subplot(2, 1, 2)
98  plt.semilogy(x, psf[:, psf.shape[0] / 2], color="blue")
99  plt.semilogy(x, psfs[:, psf.shape[0] / 2], color="red")
100  plt.xlabel("Y-axis angular distance [units of lambda/D]")
101  plt.ylabel("Normalized intensity")
102 
103  #plt.semilogy(x,psf_compass[psf.shape[0]/2,:],color="red")
104  plt.legend(["PSF exp", "PSF model"])
105  '''
106  if(correction):
107  #plt.semilogy(x,psfc[psf.shape[0]/2,:],color="purple")
108  if(synth):
109  plt.semilogy(x,psfs[psf.shape[0]/2,:],color="black")
110  #plt.legend(["PSF rec","PSF ind. assumption", "PSF COMPASS", "PSF corrected", "PSF synth"])
111  plt.legend(["PSF rec","PSF ind. assumption", "PSF COMPASS", "PSF synth"])
112  else:
113  plt.legend(["PSF rec","PSF ind. assumption", "PSF COMPASS", "PSF corrected"])
114  elif(synth):
115  plt.semilogy(x,psfs[psf.shape[0]/2,:],color="purple")
116  plt.legend(["PSF rec","PSF ind. assumption", "PSF COMPASS", "PSF synth"])
117 
118  else:
119  plt.legend(["PSF rec","PSF ind. assumption", "PSF COMPASS"])
120  '''
121 
122  f.close()
123  return psf_compass, psf, psfs
124 
125 
Here is the call graph for this function:

◆ compute_covariance_model()

def correlation_study.compute_covariance_model (   filename)

Definition at line 173 of file correlation_study.py.

173 def compute_covariance_model(filename):
174  f = h5py.File(filename, 'r')
175  Lambda_tar = f.attrs["target.Lambda"][0]
176  Lambda_wfs = f.attrs["wfs.Lambda"]
177  dt = f.attrs["ittime"]
178  gain = f.attrs["gain"]
179  wxpos = f.attrs["wfs.xpos"][0]
180  wypos = f.attrs["wfs.ypos"][0]
181  r0 = f.attrs["r0"] * (Lambda_tar / Lambda_wfs)**(6. / 5.)
182  RASC = 180. / np.pi * 3600.
183  xpos = f["dm.xpos"][:]
184  ypos = f["dm.ypos"][:]
185  p2m = f.attrs["tel_diam"] / f.attrs["pupdiam"]
186  pupshape = int(2**np.ceil(np.log2(f.attrs["pupdiam"]) + 1))
187  xactu = (xpos - pupshape / 2) * p2m
188  yactu = (ypos - pupshape / 2) * p2m
189  Ccov = np.zeros((xpos.size, xpos.size))
190  Caniso = np.zeros((xpos.size, xpos.size))
191  Cbp = np.zeros((xpos.size, xpos.size))
192  xx = np.tile(xactu, (xactu.shape[0], 1))
193  yy = np.tile(yactu, (yactu.shape[0], 1))
194  xij = xx - xx.T
195  yij = yy - yy.T
196 
197  for l in range(f.attrs["nscreens"]):
198  H = f.attrs["atm.alt"][l]
199  L0 = f.attrs["L0"][l]
200  speed = f.attrs["windspeed"][l]
201  theta = f.attrs["winddir"][l] * np.pi / 180.
202  frac = f.attrs["frac"][l]
203 
204  Htheta = np.linalg.norm([wxpos, wypos]) / RASC * H
205  vdt = speed * dt / gain
206  # Covariance matrices models on actuators space
207  M = np.zeros((xpos.size, xpos.size))
208  Mvdt = M.copy()
209  Mht = M.copy()
210  Mhvdt = M.copy()
211  angleht = np.arctan2(wypos, wxpos)
212  fc = xactu[1] - xactu[0]
213  #fc = 0.05
214 
215  M = np.linalg.norm([xij, yij], axis=0)
216  Mvdt = np.linalg.norm([xij - vdt * np.cos(theta), yij - vdt * np.sin(theta)],
217  axis=0)
218  Mht = np.linalg.norm(
219  [xij - Htheta * np.cos(angleht), yij - Htheta * np.sin(angleht)], axis=0)
220  Mhvdt = np.linalg.norm([
221  xij - vdt * np.cos(theta) - Htheta * np.cos(angleht),
222  yij - vdt * np.sin(theta) - Htheta * np.sin(angleht)
223  ], axis=0)
224  # for i in range(xpos.size):
225  # for j in range(xpos.size):
226  # Mvdt[i,j] = (np.sqrt((xactu[i]-(xactu[j]-vdt*np.cos(theta)))**2 + (yactu[i]-(yactu[j]-vdt*np.sin(theta)))**2))
227  # M[i,j] = (np.sqrt((xactu[i]-xactu[j])**2 + (yactu[i]-yactu[j])**2))
228  # Mht[i,j] = (np.sqrt((xactu[i]-(xactu[j]-Htheta*np.cos(angleht)))**2 + (yactu[i]-(yactu[j]-Htheta*np.sin(angleht)))**2))
229  # #Mhvdt[i,j] = (np.sqrt((xactu[i]-(xactu[j]+rho*np.cos(anglehvdt)))**2 + (yactu[i]-(yactu[j]+rho*np.sin(anglehvdt)))**2))
230  # Mhvdt[i,j] = (np.sqrt(((xactu[i]+vdt*np.cos(theta))-(xactu[j]-Htheta*np.cos(angleht)))**2 + ((yactu[i]+vdt*np.sin(theta))-(yactu[j]-Htheta*np.sin(angleht)))**2))
231 
232  Ccov += 0.5 * (Dphi.dphi_lowpass(Mhvdt,fc,L0,tabx,taby) - Dphi.dphi_lowpass(Mht,fc,L0,tabx,taby) \
233  - Dphi.dphi_lowpass(Mvdt,fc,L0,tabx,taby) + Dphi.dphi_lowpass(M,fc,L0,tabx,taby)) * (1./r0)**(5./3.) * frac
234 
235  Caniso += 0.5 * (Dphi.dphi_lowpass(Mht, fc, L0, tabx, taby) - Dphi.dphi_lowpass(
236  M, fc, L0, tabx, taby)) * (1. / r0)**(5. / 3.) * frac
237  Cbp += 0.5 * (Dphi.dphi_lowpass(Mvdt, fc, L0, tabx, taby) - Dphi.dphi_lowpass(
238  M, fc, L0, tabx, taby)) * (1. / r0)**(5. / 3.) * frac
239 
240  #Sp = (f.attrs["tel_diam"]/f.attrs["nxsub"])**2/2.
241  Sp = (Lambda_tar / (2 * np.pi))**2 #/3.45
242  f.close()
243  return (Caniso + Caniso.T) * Sp, (Cbp + Cbp.T) * Sp, Ccov * Sp
244 
245 
Here is the caller graph for this function:

◆ compute_psf()

def correlation_study.compute_psf (   filename)

Definition at line 26 of file correlation_study.py.

26 def compute_psf(filename):
27  otftel, otf, psf, gpu = gamora.psf_rec_Vii(filename)
28  return psf
29 
30 

◆ compute_psf_independence()

def correlation_study.compute_psf_independence (   filename)

Definition at line 31 of file correlation_study.py.

31 def compute_psf_independence(filename):
32  cov_err = rexp.get_coverr_independence(filename)
33  otfteli, otf2i, psfi, gpu = gamora.psf_rec_Vii(filenames[11], covmodes=cov_err)
34  return psfi
35 
36 

◆ cutsPSF()

def correlation_study.cutsPSF (   filename,
  psf,
  psfs 
)

Definition at line 335 of file correlation_study.py.

335 def cutsPSF(filename, psf, psfs):
336  f = h5py.File(filename, 'r')
337  Lambda_tar = f.attrs["target.Lambda"][0]
338  RASC = 180 / np.pi * 3600.
339  pixsize = Lambda_tar * 1e-6 / (
340  psf.shape[0] * f.attrs["tel_diam"] / f.attrs["pupdiam"]) * RASC
341  x = (np.arange(psf.shape[0]) - psf.shape[0] / 2) * pixsize / (
342  Lambda_tar * 1e-6 / f.attrs["tel_diam"] * RASC)
343  plt.figure()
344  plt.subplot(2, 1, 1)
345  plt.semilogy(x, psf[psf.shape[0] / 2, :], color="blue")
346  plt.semilogy(x, psfs[psf.shape[0] / 2, :], color="red")
347  plt.xlabel("X-axis angular distance [units of lambda/D]")
348  plt.ylabel("Normalized intensity")
349  plt.legend(["PSF exp", "PSF model"])
350  plt.xlim(-20, 20)
351  plt.ylim(1e-5, 1)
352  plt.subplot(2, 1, 2)
353  plt.semilogy(x, psf[:, psf.shape[0] / 2], color="blue")
354  plt.semilogy(x, psfs[:, psf.shape[0] / 2], color="red")
355  plt.xlabel("Y-axis angular distance [units of lambda/D]")
356  plt.ylabel("Normalized intensity")
357  plt.legend(["PSF exp", "PSF model"])
358  plt.xlim(-20, 20)
359  plt.ylim(1e-5, 1)
360  f.close()
361 
362 

◆ ensquare_PSF()

def correlation_study.ensquare_PSF (   filename,
  psf,
  N,
  display = False 
)

Definition at line 313 of file correlation_study.py.

313 def ensquare_PSF(filename, psf, N, display=False):
314  f = h5py.File(filename, 'r')
315  Lambda_tar = f.attrs["target.Lambda"][0]
316  RASC = 180 / np.pi * 3600.
317  pixsize = Lambda_tar * 1e-6 / (
318  psf.shape[0] * f.attrs["tel_diam"] / f.attrs["pupdiam"]) * RASC
319  x = (np.arange(psf.shape[0]) - psf.shape[0] / 2) * pixsize / (
320  Lambda_tar * 1e-6 / f.attrs["tel_diam"] * RASC)
321  w = int(N * (Lambda_tar * 1e-6 / f.attrs["tel_diam"] * RASC) / pixsize)
322  mid = psf.shape[0] / 2
323  psfe = psf[mid - w:mid + w, mid - w:mid + w]
324  if (display):
325  plt.matshow(np.log10(psfe))
326  xt = np.linspace(0, psfe.shape[0] - 1, 6).astype(np.int32)
327  yt = np.linspace(-N, N, 6).astype(np.int32)
328  plt.xticks(xt, yt)
329  plt.yticks(xt, yt)
330 
331  f.close()
332  return psf[mid - w:mid + w, mid - w:mid + w]
333 
334 

◆ filter_piston_TT()

def correlation_study.filter_piston_TT (   filename,
  C 
)

Definition at line 126 of file correlation_study.py.

126 def filter_piston_TT(filename, C):
127  IF, T = rexp.get_IF(filename)
128  IF = IF.T
129  T = T.T
130  N = IF.shape[0]
131  n = IF.shape[1]
132 
133  delta = IF.T.dot(IF).toarray() / N
134 
135  # Tip-tilt + piston
136  Tp = np.ones((T.shape[0], T.shape[1] + 1))
137  Tp[:, :2] = T.copy() #.toarray()
138  deltaT = IF.T.dot(Tp) / N
139  # Tip tilt projection on the pzt dm
140  tau = np.linalg.inv(delta).dot(deltaT)
141 
142  # Famille generatrice sans tip tilt
143  G = np.identity(n)
144  tdt = tau.T.dot(delta).dot(tau)
145  subTT = tau.dot(np.linalg.inv(tdt)).dot(tau.T).dot(delta)
146  G -= subTT
147 
148  return G.T.dot(C).dot(G)
149 
150 
Here is the caller graph for this function:

◆ filter_TT()

def correlation_study.filter_TT (   filename,
  C 
)

Definition at line 151 of file correlation_study.py.

151 def filter_TT(filename, C):
152  IF, T = rexp.get_IF(filename)
153  IF = IF.T
154  T = T.T
155  N = IF.shape[0]
156  n = IF.shape[1]
157 
158  delta = IF.T.dot(IF).toarray() / N
159 
160  deltaT = IF.T.dot(T) / N
161  # Tip tilt projection on the pzt dm
162  tau = np.linalg.inv(delta).dot(deltaT)
163 
164  # Famille generatrice sans tip tilt
165  G = np.identity(n)
166  tdt = tau.T.dot(delta).dot(tau)
167  subTT = tau.dot(np.linalg.inv(tdt)).dot(tau.T).dot(delta)
168  G -= subTT
169 
170  return G.T.dot(C).dot(G)
171 
172 
Here is the caller graph for this function:

◆ Hcor()

def correlation_study.Hcor (   f,
  Fe,
  g,
  dt 
)

Definition at line 363 of file correlation_study.py.

363 def Hcor(f, Fe, g, dt):
364  p = 1j * 2 * np.pi * f
365  return np.abs(1 / (1 + g * Fe / p * np.exp(-dt * p)))**2
366 
367 
Here is the caller graph for this function:

◆ Hretard()

def correlation_study.Hretard (   f,
  Fe,
  g,
  dt 
)

Definition at line 368 of file correlation_study.py.

368 def Hretard(f, Fe, g, dt):
369  p = 1j * 2 * np.pi * f
370  return np.abs(1 - np.exp(-p * dt / g))**2
371 
372 
Here is the caller graph for this function:

◆ load_datas()

def correlation_study.load_datas (   files)

Definition at line 263 of file correlation_study.py.

263 def load_datas(files):
264  nmodes = (files[0])["P"][:].shape[0]
265  P = (files[0])["P"][:]
266  xpos = files[0].attrs["wfs.xpos"][0]
267  ypos = files[0].attrs["wfs.ypos"][0]
268  contributors = ["tomography", "bandwidth"]
269  Lambda_tar = files[0].attrs["target.Lambda"][0]
270  Lambda_wfs = files[0].attrs["wfs.Lambda"][0]
271  L0 = files[0].attrs["L0"][0]
272  dt = files[0].attrs["ittime"]
273  H = files[0].attrs["atm.alt"][0]
274  RASC = 180 / np.pi * 3600.
275  Htheta = np.linalg.norm(
276  [xpos, ypos]
277  ) / RASC * H # np.sqrt(2)*4/RASC*H # Hardcoded for angular separation of sqrt(2)*4 arcsec
278  r0 = files[0].attrs["r0"] * (Lambda_tar / Lambda_wfs)**(6. / 5.)
279  nfiles = len(files)
280  vartomo = np.zeros((nfiles, nmodes))
281  varbp = np.zeros((nfiles, nmodes))
282  vartot = np.zeros((nfiles, nmodes))
283  theta = np.zeros(nfiles)
284  speeds = np.zeros(nfiles)
285  gain = np.zeros(nfiles)
286  # data[:,0,i] = var(tomo+bp) for file #i
287  # data[:,1,i] = var(tomo) for file #i
288  # data[:,2,i] = var(bp) for file #i
289  # data[:,3,i] = var(tomo)+var(bp) for file #i
290  ind = 0
291  print("Loading data...")
292  for f in files:
293  vartot[ind, :] = rexp.variance(f, contributors) * ((2 * np.pi / Lambda_tar)**2)
294  vartomo[ind, :] = rexp.variance(f, ["tomography"]) * (
295  (2 * np.pi / Lambda_tar)**2)
296  varbp[ind, :] = rexp.variance(f, ["bandwidth"]) * ((2 * np.pi / Lambda_tar)**2)
297  theta[ind] = f.attrs["winddir"][0]
298  speeds[ind] = f.attrs["windspeed"][0]
299  gain[ind] = float('%.1f' % f.attrs["gain"][0])
300  ind += 1
301  print(ind, "/", len(files))
302 
303  covar = (vartot - (vartomo + varbp)) / 2.
304 
305  stot = np.sum(vartot, axis=1)
306  sbp = np.sum(varbp, axis=1)
307  stomo = np.sum(vartomo, axis=1)
308  scov = np.sum(covar, axis=1)
309 
310  return stot, sbp, stomo, scov, covar
311 
312 

Variable Documentation

◆ c

correlation_study.c = colors[i]

Definition at line 487 of file correlation_study.py.

◆ color

correlation_study.color

Definition at line 490 of file correlation_study.py.

◆ colors

list correlation_study.colors = ["blue", "red", "green"]

Definition at line 480 of file correlation_study.py.

◆ datapath

string correlation_study.datapath = '/home/fferreira/Data/correlation/'

Definition at line 394 of file correlation_study.py.

◆ f

correlation_study.f = h5py.File('corStudy_Nact.h5', 'r')

Definition at line 472 of file correlation_study.py.

◆ filenames

correlation_study.filenames = glob.glob(datapath + 'roket_8m_1layer_dir*_cpu.h5')

Definition at line 395 of file correlation_study.py.

◆ files

list correlation_study.files = []

Definition at line 396 of file correlation_study.py.

◆ font

dictionary correlation_study.font = {'family': 'normal', 'weight': 'bold', 'size': 22}

Definition at line 21 of file correlation_study.py.

◆ fROKET

correlation_study.fROKET = h5py.File('ROKETStudy.h5', 'r')

Definition at line 409 of file correlation_study.py.

◆ gain

correlation_study.gain = np.zeros(nfiles)

Definition at line 405 of file correlation_study.py.

◆ ind

correlation_study.ind = 0

Definition at line 416 of file correlation_study.py.

◆ k

int correlation_study.k = 0

Definition at line 482 of file correlation_study.py.

◆ loc

correlation_study.loc

Definition at line 491 of file correlation_study.py.

◆ nfiles

correlation_study.nfiles = len(filenames)

Definition at line 402 of file correlation_study.py.

◆ nrj

correlation_study.nrj = f["nrj5"][:]

Definition at line 475 of file correlation_study.py.

◆ nrjcompass

correlation_study.nrjcompass = np.zeros(nfiles)

Definition at line 412 of file correlation_study.py.

◆ nrji

correlation_study.nrji = np.zeros(nfiles)

Definition at line 414 of file correlation_study.py.

◆ nrjroket

correlation_study.nrjroket = np.zeros(nfiles)

Definition at line 413 of file correlation_study.py.

◆ nrjs

correlation_study.nrjs = f["nrj5s"][:]

Definition at line 476 of file correlation_study.py.

◆ psf

correlation_study.psf = f["psf"][:]

Definition at line 473 of file correlation_study.py.

◆ psfi

correlation_study.psfi = fROKET["psfi"][:]

Definition at line 411 of file correlation_study.py.

◆ psfr

correlation_study.psfr = fROKET["psf"][:]

Definition at line 410 of file correlation_study.py.

◆ psfs

correlation_study.psfs = f["psfs"][:]

Definition at line 474 of file correlation_study.py.

◆ s

correlation_study.s

Definition at line 490 of file correlation_study.py.

◆ speeds

correlation_study.speeds = np.zeros(nfiles)

Definition at line 404 of file correlation_study.py.

◆ SR

correlation_study.SR = np.max(psf, axis=(0, 1))

Definition at line 477 of file correlation_study.py.

◆ SRcompass

correlation_study.SRcompass = np.zeros(nfiles)

Definition at line 406 of file correlation_study.py.

◆ SRi

correlation_study.SRi = np.zeros(nfiles)

Definition at line 408 of file correlation_study.py.

◆ SRroket

correlation_study.SRroket = np.zeros(nfiles)

Definition at line 407 of file correlation_study.py.

◆ SRs

correlation_study.SRs = np.max(psfs, axis=(0, 1))

Definition at line 478 of file correlation_study.py.

◆ tabx

correlation_study.tabx

Definition at line 400 of file correlation_study.py.

◆ taby

correlation_study.taby

Definition at line 400 of file correlation_study.py.

◆ theta

correlation_study.theta = np.zeros(nfiles)

Definition at line 403 of file correlation_study.py.

◆ v

correlation_study.v = np.unique(speeds)[i]

Definition at line 488 of file correlation_study.py.