Skip to content

Commit dafb14a

Browse files
authored
Merge pull request #107 from ajeldorado/algodiff_aj
jac: Fix Jacobian calculation wrapper to allow actuator index in the …
2 parents 3369a53 + 66bfc8b commit dafb14a

2 files changed

Lines changed: 68 additions & 49 deletions

File tree

falco/model/jacobians.py

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -121,12 +121,12 @@ def precomp(mp):
121121
pupil = pad_crop(mp.P1.compact.mask, NdmPad)
122122
Ein = pad_crop(Ein, NdmPad)
123123

124-
# Re-image the apodizer from pupil P3 back to pupil P2.
125-
if mp.flagApod:
126-
apodReimaged = pad_crop(mp.P3.compact.mask, NdmPad)
127-
apodReimaged = fp.relay(apodReimaged, NrelayFactor*mp.Nrelay2to3, mp.centering)
128-
else:
129-
apodReimaged = np.ones((NdmPad, NdmPad))
124+
# # Re-image the apodizer from pupil P3 back to pupil P2.
125+
# if mp.flagApod:
126+
# apodReimaged = pad_crop(mp.P3.compact.mask, NdmPad)
127+
# apodReimaged = fp.relay(apodReimaged, NrelayFactor*mp.Nrelay2to3, mp.centering)
128+
# else:
129+
# apodReimaged = np.ones((NdmPad, NdmPad))
130130

131131
# Compute the DM surfaces for the current DM commands
132132
if any(mp.dm_ind == 1):
@@ -190,7 +190,7 @@ def precomp(mp):
190190
NboxPad2AS = int(mp.dm2.compact.NboxAS)
191191
mp.dm2.compact.xy_box_lowerLeft_AS = mp.dm2.compact.xy_box_lowerLeft - (NboxPad2AS-mp.dm2.compact.Nbox)/2 # Account for the padding of the influence function boxes
192192

193-
apodReimaged = pad_crop(apodReimaged, mp.dm2.compact.NdmPad)
193+
# apodReimaged = pad_crop(apodReimaged, mp.dm2.compact.NdmPad)
194194
DM2stopPad = pad_crop(DM2stop, mp.dm2.compact.NdmPad)
195195
Edm2WFEpad = pad_crop(Edm2WFE, mp.dm2.compact.NdmPad)
196196

falco/model/models.py

Lines changed: 61 additions & 42 deletions
Original file line numberDiff line numberDiff line change
@@ -1124,50 +1124,62 @@ def jacobian(mp):
11241124
# for p in processes:
11251125
# p.join()
11261126

1127-
# results_Jac = [output.get() for p in processes]
1127+
# results_tuple = [output.get() for p in processes]
11281128

11291129
# # Reorder Jacobian by mode and DM from the list
11301130
# for ii in range(mp.jac.Nmode*mp.dm_ind.size):
11311131
# imode = results_order[ii][0]
11321132
# idm = results_order[ii][1]
11331133
# if idm == 1:
1134-
# jacStruct.G1[:, :, imode] = results_Jac[ii]
1134+
# jacStruct.G1[:, :, imode] = results_tuple[ii]
11351135
# if idm == 2:
1136-
# jacStruct.G2[:, :, imode] = results_Jac[ii]
1136+
# jacStruct.G2[:, :, imode] = results_tuple[ii]
11371137

11381138
print('Computing control Jacobian matrices in parallel...', end='')
11391139
# pool = multiprocessing.Pool(processes=mp.Nthreads)
11401140

11411141
with falco.util.TicToc():
1142-
results_order = [(imode, idm) for idm in mp.dm_ind for imode in range(mp.jac.Nmode)]
1143-
1144-
# # OLD WAY: with multiprocessing.Pool.starmap()
1145-
# results = pool.starmap(
1146-
# _jac_middle_layer,
1147-
# [(mp, imode, idm)for imode, idm in zip(*map(np.ravel, np.meshgrid(np.arange(mp.jac.Nmode, dtype=int), mp.dm_ind)))])
1148-
# results_Jac = results
1149-
# pool.close()
1150-
# pool.join()
1142+
# results_order = [(imode, idm) for idm in mp.dm_ind for imode in range(mp.jac.Nmode)]
1143+
1144+
# # # OLD WAY: with multiprocessing.Pool.starmap()
1145+
# # results = pool.starmap(
1146+
# # _jac_middle_layer,
1147+
# # [(mp, imode, idm)for imode, idm in zip(*map(np.ravel, np.meshgrid(np.arange(mp.jac.Nmode, dtype=int), mp.dm_ind)))])
1148+
# # results_tuple = results
1149+
# # pool.close()
1150+
# # pool.join()
1151+
1152+
# with PoolExecutor(max_workers=mp.Nthreads) as executor:
1153+
# result = executor.map(
1154+
# lambda p: _jac_middle_layer(*p),
1155+
# [(mp, imode, idm, iact)for imode, idm in zip(*map(np.ravel, np.meshgrid(np.arange(mp.jac.Nmode, dtype=int), mp.dm_ind)))]
1156+
# )
1157+
# results_tuple = tuple(result)
1158+
1159+
results_order_dm1 = [(imode, idm, iact) for iact in mp.dm1.act_ele for idm in [1, ] for imode in range(mp.jac.Nmode)]
1160+
results_order_dm2 = [(imode, idm, iact) for iact in mp.dm2.act_ele for idm in [2, ] for imode in range(mp.jac.Nmode)]
1161+
results_order = results_order_dm1 + results_order_dm2
11511162

11521163
with PoolExecutor(max_workers=mp.Nthreads) as executor:
11531164
result = executor.map(
11541165
lambda p: _jac_middle_layer(*p),
1155-
[(mp, imode, idm)for imode, idm in zip(*map(np.ravel, np.meshgrid(np.arange(mp.jac.Nmode, dtype=int), mp.dm_ind)))]
1166+
[(mp, imode, idm, iact) for imode, idm, iact in results_order]
11561167
)
1157-
results_Jac = tuple(result)
1168+
results_tuple = tuple(result)
11581169

11591170
# Reorder Jacobian by mode and DM from the list
1160-
for ii in range(mp.jac.Nmode*mp.dm_ind.size):
1171+
for ii in range(len(results_order)):
11611172
imode = results_order[ii][0]
11621173
idm = results_order[ii][1]
1174+
iact = results_order[ii][2]
11631175
if idm == 1:
1164-
jacStruct.G1[:, :, imode] = results_Jac[ii]
1176+
jacStruct.G1[:, iact, imode] = results_tuple[ii]
11651177
if idm == 2:
1166-
jacStruct.G2[:, :, imode] = results_Jac[ii]
1178+
jacStruct.G2[:, iact, imode] = results_tuple[ii]
11671179
if idm == 8:
1168-
jacStruct.G8[:, :, imode] = results_Jac[ii]
1180+
jacStruct.G8[:, iact, imode] = results_tuple[ii]
11691181
if idm == 9:
1170-
jacStruct.G9[:, :, imode] = results_Jac[ii]
1182+
jacStruct.G9[:, iact, imode] = results_tuple[ii]
11711183

11721184
print('done.')
11731185

@@ -1177,16 +1189,20 @@ def jacobian(mp):
11771189
for imode in range(mp.jac.Nmode):
11781190
if any(mp.dm_ind == 1):
11791191
print('mode%ddm%d...' % (imode, 1), end='')
1180-
jacStruct.G1[:, :, imode] = _jac_middle_layer(mp, imode, 1)
1192+
for index, iact in enumerate(mp.dm1.act_ele):
1193+
jacStruct.G1[:, index, imode] = _jac_middle_layer(mp, imode, 1, iact)
11811194
if any(mp.dm_ind == 2):
11821195
print('mode%ddm%d...' % (imode, 2), end='')
1183-
jacStruct.G2[:, :, imode] = _jac_middle_layer(mp, imode, 2)
1184-
# if any(mp.dm_ind == 8):
1185-
# print('mode%ddm%d...' % (imode, 8), end='')
1186-
# jacStruct.G8[:, :, imode] = _jac_middle_layer(mp, imode, 8)
1187-
# if any(mp.dm_ind == 9):
1188-
# print('mode%ddm%d...' % (imode, 9), end='')
1189-
# jacStruct.G9[:, :, imode] = _jac_middle_layer(mp, imode, 9)
1196+
for index, iact in enumerate(mp.dm2.act_ele):
1197+
jacStruct.G2[:, index, imode] = _jac_middle_layer(mp, imode, 2, iact)
1198+
if any(mp.dm_ind == 8):
1199+
print('mode%ddm%d...' % (imode, 8), end='')
1200+
for index, iact in enumerate(mp.dm8.act_ele):
1201+
jacStruct.G8[:, :, imode] = _jac_middle_layer(mp, imode, 8, iact)
1202+
if any(mp.dm_ind == 9):
1203+
print('mode%ddm%d...' % (imode, 9), end='')
1204+
for index, iact in enumerate(mp.dm9.act_ele):
1205+
jacStruct.G9[:, :, imode] = _jac_middle_layer(mp, imode, 9, iact)
11901206
print('done.')
11911207

11921208
# TIED ACTUATORS
@@ -1284,7 +1300,7 @@ def reassemble_jac_from_slices(mp):
12841300
return jacStruct
12851301

12861302

1287-
def _jac_middle_layer(mp, imode, idm):
1303+
def _jac_middle_layer(mp, imode, idm, iact):
12881304
"""
12891305
Select which optical layout's Jacobian model to use and get E-field.
12901306
@@ -1301,34 +1317,37 @@ def _jac_middle_layer(mp, imode, idm):
13011317
"""
13021318
if mp.layout.lower() in ('fourier', 'proper', 'fpm_scale'):
13031319

1304-
if idm == 1:
1305-
Nele = mp.dm1.Nele
1306-
act_ele = mp.dm1.act_ele
1307-
elif idm == 2:
1308-
Nele = mp.dm2.Nele
1309-
act_ele = mp.dm2.act_ele
1310-
jacMode = np.zeros((mp.Fend.corr.Npix, Nele), dtype=complex)
1320+
# if idm == 1:
1321+
# Nele = mp.dm1.Nele
1322+
# act_ele = mp.dm1.act_ele
1323+
# elif idm == 2:
1324+
# Nele = mp.dm2.Nele
1325+
# act_ele = mp.dm2.act_ele
1326+
# jacMode = np.zeros((mp.Fend.corr.Npix, Nele), dtype=complex)
13111327

13121328
if mp.coro.upper() in ('LC', 'APLC', 'HLC', 'FLC', 'SPLC'):
1313-
for index, iact in enumerate(act_ele):
1314-
jacMode[:, index] = jacobians.lyot(mp, imode, idm, iact)
1329+
jacModeAct = jacobians.lyot(mp, imode, idm, iact)
1330+
# for index, iact in enumerate(act_ele):
1331+
# jacMode[:, index] = jacobians.lyot(mp, imode, idm, iact)
13151332

13161333
elif mp.coro.upper() in ('VC', 'AVC', 'VORTEX'):
1317-
for index, iact in enumerate(act_ele):
1318-
jacMode[:, index] = jacobians.vortex(mp, imode, idm, iact)
1334+
jacModeAct = jacobians.vortex(mp, imode, idm, iact)
1335+
# for index, iact in enumerate(act_ele):
1336+
# jacMode[:, index] = jacobians.vortex(mp, imode, idm, iact)
13191337

13201338
elif mp.layout.lower() in ('wfirst_phaseb_proper', 'roman_phasec_proper'):
13211339
if mp.coro.upper() in ('HLC', 'SPC', 'SPLC'):
1322-
for index, iact in enumerate(act_ele):
1323-
jacMode[:, index] = jacobians.lyot(mp, imode, idm, iact)
1340+
jacModeAct = jacobians.lyot(mp, imode, idm, iact)
1341+
# for index, iact in enumerate(act_ele):
1342+
# jacMode[:, index] = jacobians.lyot(mp, imode, idm, iact)
13241343
else:
13251344
raise ValueError('%s not recognized as value for mp.coro' %
13261345
mp.coro)
13271346

13281347
else:
13291348
raise ValueError('mp.layout.lower not recognized')
13301349

1331-
return jacMode
1350+
return jacModeAct
13321351

13331352

13341353

0 commit comments

Comments
 (0)