4 sys.path.insert(0, os.environ[
"SHESHA_ROOT"] +
"/widgets/")
5 sys.path.insert(0, os.environ[
"SHESHA_ROOT"] +
"/lib/")
10 import carmaWrap
as ch
13 import matplotlib.pyplot
as plt
14 import hdf5_util
as h5u
15 import resDataBase
as db
16 import astropy.io.fits
as pf
20 print(
"TEST SHESHA\n closed loop: call loop(int niter)")
22 pathResults =
"/volumes/hra/micado/PYR39m_RoundPupil_RUN2/"
23 dBResult =
"/volumes/hra/micado/PYR_39m_RoundPupil_RUN2.h5"
26 if (len(sys.argv) == 1):
27 error =
'command line should be:"python -i test.py parameters_filename"\n with "parameters_filename" the path to the parameters file'
28 raise Exception(error)
29 if (len(sys.argv) == 2):
30 print(
"Using Internal parameters...")
44 print(
"DETECTED BASH SCRIPT")
47 freqs = [float(sys.argv[2])]
48 RONS = [float(sys.argv[3])]
49 MODU = [float(sys.argv[4])]
50 gainslist = [float(sys.argv[5])]
51 magnitudes = [float(sys.argv[6])]
52 nKL_Filt = float(sys.argv[7])
55 Nsimutot = len(gainslist) * len(freqs) * len(RONS) * len(MODU) * len(magnitudes)
57 if (
not glob.glob(pathResults)):
58 print(
"Results folder not found. Creating it now:")
59 tools.system(
"mkdir " + pathResults)
60 if (
not glob.glob(pathResults +
"PSFs/")):
61 print(
"PSFs folder not found. Creating it now:")
63 tools.system(
"mkdir " + pathResults +
"PSFs/")
66 param_file = sys.argv[1]
67 if (param_file.split(
'.')[-1] == b
"py"):
68 filename = param_file.split(
'/')[-1]
69 param_path = param_file.split(filename)[0]
70 sys.path.insert(0, param_path)
71 exec(
"import %s as config" % filename.split(
".py")[0])
73 elif (param_file.split(
'.')[-1] == b
"h5"):
74 sys.path.insert(0, os.environ[
"SHESHA_ROOT"] +
"/data/par/par4bench/")
75 import scao_sh_16x16_8pix
as config
77 h5u.configFromH5(param_file, config)
79 raise ValueError(
"Parameter file extension must be .py or .h5")
81 print(
"param_file is", param_file)
83 if (hasattr(config,
"simul_name")):
84 if (config.simul_name
is None):
87 simul_name = config.simul_name
90 print(
"simul name is", simul_name)
93 if (simul_name == b
""):
97 param_dict = h5u.params_dictionary(config)
98 matricesToLoad = h5u.checkMatricesDataBase(os.environ[
"SHESHA_ROOT"] +
"/data/",
102 c = ch.carmaWrap_context(devices=np.array([0, 1, 2, 3], dtype=np.int32))
108 hdulist = pf.open(filepath)
109 header = hdulist[0].header
110 names = np.sort(list(set(df))).tolist()
113 if (type(val)
is list):
116 value += (str(v) +
" ")
117 elif (type(val)
is np.ndarray):
120 value += (str(v) +
" ")
123 header.set(name, value,
'')
124 hdulist.writeto(filepath, clobber=
True)
129 param_dict = h5u.params_dictionary(config)
130 matricesToLoad = h5u.checkMatricesDataBase(os.environ[
"SHESHA_ROOT"] +
"/data/",
133 wfs, tel = ao.wfs_init(config.p_wfss, config.p_atmos, config.p_tel, config.p_geom,
134 config.p_target, config.p_loop, config.p_dms)
137 atm = ao.atmos_init(c, config.p_atmos, config.p_tel, config.p_geom, config.p_loop,
138 config.p_wfss, wfs, config.p_target, rank=0, clean=clean,
143 dms = ao.dm_init(config.p_dms, config.p_wfss, wfs, config.p_geom, config.p_tel)
147 tar = ao.target_init(c, tel, config.p_target, config.p_atmos, config.p_geom,
148 config.p_tel, config.p_dms)
152 rtc = ao.rtc_init(tel, wfs, config.p_wfss, dms, config.p_dms, config.p_geom,
153 config.p_rtc, config.p_atmos, atm, config.p_tel, config.p_loop,
154 do_refslp=
False, clean=clean, simul_name=simul_name,
155 load=matricesToLoad, doimat=0)
157 h5u.validDataBase(os.environ[
"SHESHA_ROOT"] +
"/data/", matricesToLoad)
159 print(
"====================")
161 print(
"====================")
162 print(
"objects initialzed on GPU:")
163 print(
"--------------------------------------------------------")
169 return wfs, tel, atm, dms, tar, rtc
172 def loop(n, wfs, tel, atm, dms, tar, rtc):
174 print(
"----------------------------------------------------")
175 print(
"iter# | S.E. SR | L.E. SR | Est. Rem. | framerate")
176 print(
"----------------------------------------------------")
183 if (config.p_controllers[0].type == b
"geo"):
184 for t
in range(config.p_target.ntargets):
185 tar.atmos_trace(t, atm, tel)
186 rtc.do_control_geo(0, dms, tar, 0)
190 for t
in range(config.p_target.ntargets):
191 tar.atmos_trace(t, atm, tel)
193 for w
in range(len(config.p_wfss)):
194 wfs.raytrace(w,
"all", tel, atm, dms)
195 wfs.sensors_compimg(w)
201 if ((i + 1) % 100 == 0):
202 print(
"Iter#:", i + 1)
205 SR = tar.get_strehl(t)
206 print(
"Tar %d at %3.2fMicrons:" % (t + 1, tar.Lambda[t]))
207 signal_se =
"SR S.E: %1.2f " % SR[0]
208 signal_le =
"SR L.E: %1.2f " % SR[1]
210 print(signal_se + signal_le)
214 numiter.append(i + 1)
225 print(
" loop execution time:", t1 - t0,
" (", n,
"iterations), ", (t1 - t0) / n,
226 "(mean) ", n / (t1 - t0),
"Hz")
228 for t
in range(config.p_target.ntargets):
229 SR = tar.get_strehl(t)
231 return SRList, tar.Lambda.tolist(), sr_le, sr_se, numiter
237 dictProcess, dictplot = getDataFrameColumns()
238 colnames = h5u.params_dictionary(config)
239 resAll = pd.DataFrame( columns=colnames.keys()) # res is the local dataframe for THIS data set
240 resAll = resAll.append(colnames, ignore_index=True) #Fill dataframe
244 colnames = h5u.params_dictionary(config)
247 "PSFFilenames":
None,
249 "lambdaTarget":
None,
264 resAll = db.readDataBase(fullpath=dBResult)
265 if (
not (type(resAll) == pd.core.frame.DataFrame)):
266 print(
"Creating compass database")
267 resAll = db.createDf(list(colnames.keys()) + list(
286 config.p_loop.set_ittime(1 / freq)
288 config.p_wfs0.set_noise(RON)
289 for modulation
in MODU:
291 config.p_wfs0.set_pyr_npts(
292 int(np.ceil(int(rMod * 2 * 3.141592653589793) / 4.) * 4))
293 config.p_wfs0.set_pyr_ampl(rMod)
294 for gain
in gainslist:
295 config.p_controller0.set_gain(gain)
296 for magnitude
in magnitudes:
298 config.p_wfs0.set_gsmag(magnitude)
300 columns=list(colnames.keys()) + list(
302 print(
"Freq = %3.2f Hz" % (1. / config.p_loop.ittime))
303 print(
"Magnitude = %3.2f" % config.p_wfs0.gsmag)
304 print(
"Gain = %3.2f" % config.p_controller0.gain)
309 print(
"Reading cMat")
313 print(
"Setting cMat")
315 rtc.set_cmat(0, cmat.copy())
316 print(
"Starting Loop")
317 SR, lambdaTargetList, sr_le, sr_se, numiter =
loop(
318 config.p_loop.niter, wfs, tel, atm, dms, tar, rtc)
319 dfparams = h5u.params_dictionary(
321 dfparams.update(simunames)
323 res = db.fillDf(res, dfparams)
324 res.loc[0,
"NklFilt"] = nKL_Filt
325 res.loc[0,
"Nkl"] = cmat.shape[0] - 2 - nKL_Filt
326 res.loc[0,
"NklTot"] = cmat.shape[0] - 2
327 res.loc[0,
"Nactu"] = cmat.shape[0]
328 res.loc[0,
"Nslopes"] = cmat.shape[1]
329 res.loc[0,
"Nphotons"] = config.p_wfs0._nphotons
330 res.loc[0,
"RON"] = RON
331 res.loc[0,
"Nphotons"] = config.p_wfs0._nphotons
333 res.srir.values[0] = SR
334 res.lambdaTarget.values[0] = lambdaTargetList
335 res.loc[0,
"gsmag"] = config.p_wfs0.gsmag
336 res.loc[0,
"gain"] = config.p_controller0.gain
340 res.loc[0,
"simulname"] = simulName
341 print(
"Saving PSFs...")
343 for t
in range(config.p_target.ntargets):
344 PSFtarget = tar.get_image(t,
"le")
345 date = time.strftime(
"_%d-%m-%Y_%H:%M:%S_")
346 lam =
"%3.2f" % tar.Lambda.tolist()[t]
347 lam = lam.replace(
".",
"_")
348 PSFName =
"PYR_" + lam +
"_" + date +
".fits"
349 PSFNameList.append(PSFName)
351 pf.writeto(pathResults +
"PSFs/" + PSFName, PSFtarget.copy(),
353 lam2 =
"%3.2f" % tar.Lambda.tolist()[t]
354 res.loc[0,
"SR_%s" % lam2] = SR[t]
355 filepath = pathResults +
"PSFs/" + PSFName
358 hdulist = pf.open(filepath)
359 header = hdulist[0].header
361 header[
"wavelength"] = tar.Lambda.tolist()[t]
362 hdulist.writeto(filepath, clobber=
True)
366 res.loc[0,
"type_ap"] = str(res.loc[0,
"type_ap"][0])
367 res.loc[0,
"type"] = str(res.loc[0,
"type"][0])
368 res.loc[0,
"type"] =
"pzt, tt"
369 res.loc[0,
"npix"] = res.loc[0,
"npix"][0]
371 res.loc[0,
"pixsize"] = res.loc[0,
"pixsize"][0]
372 res.PSFFilenames.values[0] = PSFNameList
373 resAll = db.fillDf(resAll, res)
375 resAll.to_hdf(dBResult,
"resAll", complevel=9, complib=
'blosc')
378 print(
"Simulation Done...")
380 Sauver PSF dans le bon nom + directory
381 ranger... + params dans le header
383 SR = np.zeros((3, len(set(resAll.gsmag))))
385 for mag in list(set(resAll.gsmag)):
386 SR[0,i] = resAll[resAll.gsmag == mag]["SR_1.20"].max()
387 SR[1,i] = resAll[resAll.gsmag == mag]["SR_1.65"].max()
388 SR[2,i] = resAll[resAll.gsmag == mag]["SR_2.20"].max()