Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1#!/usr/bin/env python 

2# encoding: utf-8 

3""" 

4*detect arc-lines on a pinhole frame to generate a dispersion solution* 

5 

6:Author: 

7 Marco Landoni & David Young 

8 

9:Date Created: 

10 September 1, 2020 

11""" 

12################# GLOBAL IMPORTS #################### 

13from builtins import object 

14import sys 

15import os 

16os.environ['TERM'] = 'vt100' 

17from fundamentals import tools 

18from soxspipe.commonutils import keyword_lookup 

19from soxspipe.commonutils import detector_lookup 

20from astropy.stats import sigma_clipped_stats 

21from photutils import datasets 

22from photutils import DAOStarFinder 

23from scipy.optimize import curve_fit 

24from fundamentals.renderer import list_of_dictionaries 

25from os.path import expanduser 

26import matplotlib.pyplot as plt 

27from os.path import expanduser 

28from astropy.stats import sigma_clip, mad_std 

29import numpy as np 

30import math 

31from photutils.utils import NoDetectionsWarning 

32import warnings 

33from astropy.visualization import hist 

34from soxspipe.commonutils.polynomials import chebyshev_order_wavelength_polynomials 

35from soxspipe.commonutils.filenamer import filenamer 

36from soxspipe.commonutils.dispersion_map_to_pixel_arrays import dispersion_map_to_pixel_arrays 

37import pandas as pd 

38from tabulate import tabulate 

39 

40 

41class create_dispersion_map(object): 

42 """ 

43 *detect arc-lines on a pinhole frame to generate a dispersion solution* 

44 

45 **Key Arguments:** 

46 - ``log`` -- logger 

47 - ``settings`` -- the settings dictionary 

48 - ``pinholeFrame`` -- the calibrated pinhole frame (single or multi) 

49 - ``firstGuessMap`` -- the first guess dispersion map from the `soxs_disp_solution` recipe (needed in `soxs_spat_solution` recipe). Default *False*. 

50 - ``qcTable`` -- the data frame to collect measured QC metrics 

51 - ``productsTable`` -- the data frame to collect output products 

52 

53 **Usage:** 

54 

55 ```python 

56 from soxspipe.commonutils import create_dispersion_map 

57 mapPath = create_dispersion_map( 

58 log=log, 

59 settings=settings, 

60 pinholeFrame=frame, 

61 firstGuessMap=False 

62 ).get() 

63 ``` 

64 """ 

65 

66 def __init__( 

67 self, 

68 log, 

69 settings, 

70 pinholeFrame, 

71 firstGuessMap=False, 

72 qcTable=False, 

73 productsTable=False 

74 ): 

75 self.log = log 

76 log.debug("instantiating a new 'create_dispersion_map' object") 

77 self.settings = settings 

78 self.pinholeFrame = pinholeFrame 

79 self.firstGuessMap = firstGuessMap 

80 self.qc = qcTable 

81 self.products = productsTable 

82 

83 # KEYWORD LOOKUP OBJECT - LOOKUP KEYWORD FROM DICTIONARY IN RESOURCES 

84 # FOLDER 

85 kw = keyword_lookup( 

86 log=self.log, 

87 settings=self.settings 

88 ).get 

89 self.kw = kw 

90 self.arm = pinholeFrame.header[kw("SEQ_ARM")] 

91 self.dateObs = pinholeFrame.header[kw("DATE_OBS")] 

92 

93 # DETECTOR PARAMETERS LOOKUP OBJECT 

94 self.detectorParams = detector_lookup( 

95 log=log, 

96 settings=settings 

97 ).get(self.arm) 

98 

99 warnings.simplefilter('ignore', NoDetectionsWarning) 

100 

101 return None 

102 

103 def get(self): 

104 """ 

105 *generate the dispersion map* 

106 

107 **Return:** 

108 - ``mapPath`` -- path to the file containing the coefficients of the x,y polynomials of the global dispersion map fit 

109 """ 

110 self.log.debug('starting the ``get`` method') 

111 

112 # WHICH RECIPE ARE WE WORKING WITH? 

113 if self.firstGuessMap: 

114 recipe = "soxs-spatial-solution" 

115 slit_deg = self.settings[recipe]["slit-deg"] 

116 else: 

117 recipe = "soxs-disp-solution" 

118 slit_deg = 0 

119 

120 # READ PREDICTED LINE POSITIONS FROM FILE - RETURNED AS DATAFRAME 

121 orderPixelTable = self.get_predicted_line_list() 

122 

123 # GET THE WINDOW SIZE FOR ATTEMPTING TO DETECT LINES ON FRAME 

124 windowSize = self.settings[recipe]["pixel-window-size"] 

125 self.windowHalf = int(windowSize / 2) 

126 

127 # DETECT THE LINES ON THE PINHILE FRAME AND 

128 # ADD OBSERVED LINES TO DATAFRAME 

129 orderPixelTable = orderPixelTable.apply( 

130 self.detect_pinhole_arc_line, axis=1) 

131 

132 # DROP MISSING VALUES 

133 orderPixelTable.dropna(axis='index', how='any', subset=[ 

134 'observed_x'], inplace=True) 

135 

136 order_deg = self.settings[recipe]["order-deg"] 

137 wavelength_deg = self.settings[ 

138 recipe]["wavelength-deg"] 

139 

140 # ITERATIVELY FIT THE POLYNOMIAL SOLUTIONS TO THE DATA 

141 popt_x, popt_y = self.fit_polynomials( 

142 orderPixelTable=orderPixelTable, 

143 wavelength_deg=wavelength_deg, 

144 order_deg=order_deg, 

145 slit_deg=slit_deg, 

146 ) 

147 

148 # WRITE THE MAP TO FILE 

149 mapPath = self.write_map_to_file( 

150 popt_x, popt_y, order_deg, wavelength_deg, slit_deg) 

151 

152 self.log.debug('completed the ``get`` method') 

153 return mapPath 

154 

155 def get_predicted_line_list( 

156 self): 

157 """*lift the predicted line list from the static calibrations* 

158 

159 **Return:** 

160 - ``orderPixelTable`` -- a panda's data-frame containing wavelength,order,slit_index,slit_position,detector_x,detector_y 

161 """ 

162 self.log.debug('starting the ``get_predicted_line_list`` method') 

163 

164 kw = self.kw 

165 pinholeFrame = self.pinholeFrame 

166 dp = self.detectorParams 

167 

168 # WHICH TYPE OF PINHOLE FRAME DO WE HAVE - SINGLE OR MULTI 

169 if self.pinholeFrame.header[kw("DPR_TECH")] == "ECHELLE,PINHOLE": 

170 frameTech = "single" 

171 elif self.pinholeFrame.header[kw("DPR_TECH")] == "ECHELLE,MULTI-PINHOLE": 

172 frameTech = "multi" 

173 else: 

174 raise TypeError( 

175 "The input frame needs to be a calibrated single- or multi-pinhole arc lamp frame") 

176 

177 # FIND THE APPROPRIATE PREDICTED LINE-LIST 

178 arm = self.arm 

179 if kw('WIN_BINX') in pinholeFrame.header: 

180 binx = int(self.pinholeFrame.header[kw('WIN_BINX')]) 

181 biny = int(self.pinholeFrame.header[kw('WIN_BINY')]) 

182 else: 

183 binx = 1 

184 biny = 1 

185 

186 # READ THE FILE 

187 home = expanduser("~") 

188 calibrationRootPath = self.settings[ 

189 "calibration-data-root"].replace("~", home) 

190 predictedLinesFile = calibrationRootPath + "/" + dp["predicted pinhole lines"][frameTech][f"{binx}x{biny}"] 

191 

192 # READ CSV FILE TO PANDAS DATAFRAME 

193 orderPixelTable = pd.read_csv(predictedLinesFile) 

194 

195 # RENAME ALL COLUMNS FOR CONSISTENCY 

196 listName = [] 

197 listName[:] = [l if l else l for l in listName] 

198 orderPixelTable.columns = [d.lower() if d.lower() in [ 

199 "order", "wavelength"] else d for d in orderPixelTable.columns] 

200 

201 # WANT TO DETERMINE SYSTEMATIC SHIFT IF FIRST GUESS SOLUTION PRESENT 

202 if self.firstGuessMap: 

203 # ADD SOME EXTRA COLUMNS TO DATAFRAME 

204 

205 # FILTER THE PREDICTED LINES TO ONLY SLIT POSITION INCLUDED IN 

206 # SINGLE PINHOLE FRAMES 

207 slitIndex = int(dp["mid_slit_index"]) 

208 

209 # GET THE OBSERVED PIXELS VALUES 

210 orderPixelTable = dispersion_map_to_pixel_arrays( 

211 log=self.log, 

212 dispersionMapPath=self.firstGuessMap, 

213 orderPixelTable=orderPixelTable 

214 ) 

215 

216 # CREATE A COPY OF THE DATA-FRAME TO DETERMINE SHIFTS 

217 tmpList = orderPixelTable.copy() 

218 

219 mask = (tmpList['slit_index'] == slitIndex) 

220 tmpList.loc[mask, 'shift_x'] = tmpList.loc[ 

221 mask, 'detector_x'].values - tmpList.loc[mask, 'fit_x'].values 

222 tmpList.loc[mask, 'shift_y'] = tmpList.loc[ 

223 mask, 'detector_y'].values - tmpList.loc[mask, 'fit_y'].values 

224 

225 # MERGING SHIFTS INTO MAIN DATAFRAME 

226 tmpList = tmpList.loc[tmpList['shift_x'].notnull( 

227 ), ['wavelength', 'order', 'shift_x', 'shift_y']] 

228 orderPixelTable = orderPixelTable.merge(tmpList, on=[ 

229 'wavelength', 'order'], how='outer') 

230 

231 # DROP ROWS WITH MISSING SHIFTS 

232 orderPixelTable.dropna(axis='index', how='any', subset=[ 

233 'shift_x'], inplace=True) 

234 

235 # SHIFT DETECTOR LINE PIXEL POSITIONS BY SHIFTS 

236 # UPDATE FILTERED VALUES 

237 orderPixelTable.loc[ 

238 :, 'detector_x'] -= orderPixelTable.loc[:, 'shift_x'] 

239 orderPixelTable.loc[ 

240 :, 'detector_y'] -= orderPixelTable.loc[:, 'shift_y'] 

241 

242 # DROP HELPER COLUMNS 

243 orderPixelTable.drop(columns=['fit_x', 'fit_y', 

244 'shift_x', 'shift_y'], inplace=True) 

245 

246 self.log.debug('completed the ``get_predicted_line_list`` method') 

247 return orderPixelTable 

248 

249 def detect_pinhole_arc_line( 

250 self, 

251 predictedLine): 

252 """*detect the observed position of an arc-line given the predicted pixel positions* 

253 

254 **Key Arguments:** 

255 - ``predictedLine`` -- single predicted line coordinates from predicted line-list 

256 

257 **Return:** 

258 - ``predictedLine`` -- the line with the observed pixel coordinates appended (if detected, otherwise nan) 

259 """ 

260 self.log.debug('starting the ``detect_pinhole_arc_line`` method') 

261 

262 pinholeFrame = self.pinholeFrame 

263 windowHalf = self.windowHalf 

264 x = predictedLine['detector_x'] 

265 y = predictedLine['detector_y'] 

266 

267 # CLIP A STAMP FROM IMAGE AROUNDS PREDICTED POSITION 

268 xlow = int(np.max([x - windowHalf, 0])) 

269 xup = int(np.min([x + windowHalf, pinholeFrame.shape[1]])) 

270 ylow = int(np.max([y - windowHalf, 0])) 

271 yup = int(np.min([y + windowHalf, pinholeFrame.shape[0]])) 

272 stamp = pinholeFrame[ylow:yup, xlow:xup] 

273 # CONVERT TO MASKED ARRAY 

274 stamp = np.asanyarray(stamp) 

275 

276 # USE DAOStarFinder TO FIND LINES WITH 2D GUASSIAN FITTING 

277 mean, median, std = sigma_clipped_stats(stamp, sigma=3.0) 

278 

279 daofind = DAOStarFinder( 

280 fwhm=2.0, threshold=5. * std, roundlo=-3.0, roundhi=3.0, sharplo=-3.0, sharphi=3.0) 

281 sources = daofind(stamp - median) 

282 

283 # plt.clf() 

284 # plt.imshow(stamp) 

285 old_resid = windowHalf * 4 

286 if sources: 

287 # FIND SOURCE CLOSEST TO CENTRE 

288 if len(sources) > 1: 

289 for source in sources: 

290 tmp_x = source['xcentroid'] 

291 tmp_y = source['ycentroid'] 

292 new_resid = ((windowHalf - tmp_x)**2 + 

293 (windowHalf - tmp_y)**2)**0.5 

294 if new_resid < old_resid: 

295 observed_x = tmp_x + xlow 

296 observed_y = tmp_y + ylow 

297 old_resid = new_resid 

298 else: 

299 observed_x = sources[0]['xcentroid'] + xlow 

300 observed_y = sources[0]['ycentroid'] + ylow 

301 # plt.scatter(observed_x - xlow, observed_y - 

302 # ylow, marker='x', s=30) 

303 # plt.show() 

304 else: 

305 observed_x = np.nan 

306 observed_y = np.nan 

307 # plt.show() 

308 

309 predictedLine['observed_x'] = observed_x 

310 predictedLine['observed_y'] = observed_y 

311 

312 self.log.debug('completed the ``detect_pinhole_arc_line`` method') 

313 return predictedLine 

314 

315 def write_map_to_file( 

316 self, 

317 xcoeff, 

318 ycoeff, 

319 order_deg, 

320 wavelength_deg, 

321 slit_deg): 

322 """*write out the fitted polynomial solution coefficients to file* 

323 

324 **Key Arguments:** 

325 - ``xcoeff`` -- the x-coefficients 

326 - ``ycoeff`` -- the y-coefficients 

327 - ``order_deg`` -- degree of the order fitting 

328 - ``wavelength_deg`` -- degree of wavelength fitting 

329 - ``slit_deg`` -- degree of the slit fitting (False for single pinhole) 

330 

331 **Return:** 

332 - ``disp_map_path`` -- path to the saved file 

333 """ 

334 self.log.debug('starting the ``write_map_to_file`` method') 

335 

336 arm = self.arm 

337 

338 # SORT X COEFFICIENT OUTPUT TO WRITE TO FILE 

339 coeff_dict_x = {} 

340 coeff_dict_x["axis"] = "x" 

341 coeff_dict_x["order-deg"] = order_deg 

342 coeff_dict_x["wavelength-deg"] = wavelength_deg 

343 coeff_dict_x["slit-deg"] = slit_deg 

344 n_coeff = 0 

345 for i in range(0, order_deg + 1): 

346 for j in range(0, wavelength_deg + 1): 

347 for k in range(0, slit_deg + 1): 

348 coeff_dict_x[f'c{i}{j}{k}'] = xcoeff[n_coeff] 

349 n_coeff += 1 

350 

351 # SORT Y COEFFICIENT OUTPUT TO WRITE TO FILE 

352 coeff_dict_y = {} 

353 coeff_dict_y["axis"] = "y" 

354 coeff_dict_y["order-deg"] = order_deg 

355 coeff_dict_y["wavelength-deg"] = wavelength_deg 

356 coeff_dict_y["slit-deg"] = slit_deg 

357 n_coeff = 0 

358 for i in range(0, order_deg + 1): 

359 for j in range(0, wavelength_deg + 1): 

360 for k in range(0, slit_deg + 1): 

361 coeff_dict_y[f'c{i}{j}{k}'] = ycoeff[n_coeff] 

362 n_coeff += 1 

363 

364 # DETERMINE WHERE TO WRITE THE FILE 

365 home = expanduser("~") 

366 outDir = self.settings["intermediate-data-root"].replace("~", home) 

367 

368 filename = filenamer( 

369 log=self.log, 

370 frame=self.pinholeFrame, 

371 settings=self.settings 

372 ) 

373 if slit_deg == 0: 

374 filename = filename.split("ARC")[0] + "DISP_MAP.csv" 

375 else: 

376 filename = filename.split("ARC")[0] + "2D_MAP.csv" 

377 filePath = f"{outDir}/{filename}" 

378 dataSet = list_of_dictionaries( 

379 log=self.log, 

380 listOfDictionaries=[coeff_dict_x, coeff_dict_y] 

381 ) 

382 csvData = dataSet.csv(filepath=filePath) 

383 

384 self.log.debug('completed the ``write_map_to_file`` method') 

385 return filePath 

386 

387 def calculate_residuals( 

388 self, 

389 orderPixelTable, 

390 xcoeff, 

391 ycoeff, 

392 order_deg, 

393 wavelength_deg, 

394 slit_deg): 

395 """*calculate residuals of the polynomial fits against the observed line positions* 

396 

397 **Key Arguments:** 

398 

399 - ``orderPixelTable`` -- the predicted line list as a data frame 

400 - ``xcoeff`` -- the x-coefficients 

401 - ``ycoeff`` -- the y-coefficients 

402 - ``order_deg`` -- degree of the order fitting 

403 - ``wavelength_deg`` -- degree of wavelength fitting 

404 - ``slit_deg`` -- degree of the slit fitting (False for single pinhole) 

405 

406 **Return:** 

407 - ``residuals`` -- combined x-y residuals 

408 - ``mean`` -- the mean of the combine residuals 

409 - ``std`` -- the stdev of the combine residuals 

410 - ``median`` -- the median of the combine residuals 

411 """ 

412 self.log.debug('starting the ``calculate_residuals`` method') 

413 

414 arm = self.arm 

415 

416 # POLY FUNCTION NEEDS A DATAFRAME AS INPUT 

417 poly = chebyshev_order_wavelength_polynomials( 

418 log=self.log, order_deg=order_deg, wavelength_deg=wavelength_deg, slit_deg=slit_deg).poly 

419 

420 # CALCULATE X & Y RESIDUALS BETWEEN OBSERVED LINE POSITIONS AND POLY 

421 # FITTED POSITIONS 

422 orderPixelTable["fit_x"] = poly(orderPixelTable, *xcoeff) 

423 orderPixelTable["fit_y"] = poly(orderPixelTable, *ycoeff) 

424 orderPixelTable["residuals_x"] = orderPixelTable[ 

425 "fit_x"] - orderPixelTable["observed_x"] 

426 orderPixelTable["residuals_y"] = orderPixelTable[ 

427 "fit_y"] - orderPixelTable["observed_y"] 

428 

429 # CALCULATE COMBINED RESIDUALS AND STATS 

430 orderPixelTable["residuals_xy"] = np.sqrt(np.square( 

431 orderPixelTable["residuals_x"]) + np.square(orderPixelTable["residuals_y"])) 

432 combined_res_mean = np.mean(orderPixelTable["residuals_xy"]) 

433 combined_res_std = np.std(orderPixelTable["residuals_xy"]) 

434 combined_res_median = np.median(orderPixelTable["residuals_xy"]) 

435 

436 self.log.debug('completed the ``calculate_residuals`` method') 

437 return combined_res_mean, combined_res_std, combined_res_median, orderPixelTable 

438 

439 def fit_polynomials( 

440 self, 

441 orderPixelTable, 

442 wavelength_deg, 

443 order_deg, 

444 slit_deg): 

445 """*iteratively fit the dispersion map polynomials to the data, clipping residuals with each iteration* 

446 

447 **Key Arguments:** 

448 - ``orderPixelTable`` -- data frame containing order, wavelengths, slit positions and observed pixel positions 

449 - ``wavelength_deg`` -- degree of wavelength fitting 

450 - ``order_deg`` -- degree of the order fitting 

451 - ``slit_deg`` -- degree of the slit fitting (0 for single pinhole) 

452 

453 **Return:** 

454 - ``xcoeff`` -- the x-coefficients post clipping 

455 - ``ycoeff`` -- the y-coefficients post clipping 

456 """ 

457 self.log.debug('starting the ``fit_polynomials`` method') 

458 

459 arm = self.arm 

460 

461 if self.firstGuessMap: 

462 recipe = "soxs-spatial-solution" 

463 else: 

464 recipe = "soxs-disp-solution" 

465 

466 clippedCount = 1 

467 

468 poly = chebyshev_order_wavelength_polynomials( 

469 log=self.log, order_deg=order_deg, wavelength_deg=wavelength_deg, slit_deg=slit_deg).poly 

470 

471 clippingSigma = self.settings[ 

472 recipe]["poly-fitting-residual-clipping-sigma"] 

473 clippingIterationLimit = self.settings[ 

474 recipe]["clipping-iteration-limit"] 

475 

476 iteration = 0 

477 while clippedCount > 0 and iteration < clippingIterationLimit: 

478 iteration += 1 

479 observed_x = orderPixelTable["observed_x"].to_numpy() 

480 observed_y = orderPixelTable["observed_y"].to_numpy() 

481 # USE LEAST-SQUARED CURVE FIT TO FIT CHEBY POLYS 

482 # FIRST X 

483 coeff = np.ones((order_deg + 1) * 

484 (wavelength_deg + 1) * (slit_deg + 1)) 

485 self.log.info("""curvefit x""" % locals()) 

486 

487 xcoeff, pcov_x = curve_fit( 

488 poly, xdata=orderPixelTable, ydata=observed_x, p0=coeff) 

489 

490 # NOW Y 

491 self.log.info("""curvefit y""" % locals()) 

492 ycoeff, pcov_y = curve_fit( 

493 poly, xdata=orderPixelTable, ydata=observed_y, p0=coeff) 

494 

495 self.log.info("""calculate_residuals""" % locals()) 

496 mean_res, std_res, median_res, orderPixelTable = self.calculate_residuals( 

497 orderPixelTable=orderPixelTable, 

498 xcoeff=xcoeff, 

499 ycoeff=ycoeff, 

500 order_deg=order_deg, 

501 wavelength_deg=wavelength_deg, 

502 slit_deg=slit_deg) 

503 

504 # SIGMA-CLIP THE DATA 

505 self.log.info("""sigma_clip""" % locals()) 

506 masked_residuals = sigma_clip( 

507 orderPixelTable["residuals_xy"], sigma_lower=clippingSigma, sigma_upper=clippingSigma, maxiters=1, cenfunc='median', stdfunc=mad_std) 

508 orderPixelTable["residuals_masked"] = masked_residuals.mask 

509 # RETURN BREAKDOWN OF COLUMN VALUE COUNT 

510 valCounts = orderPixelTable[ 

511 'residuals_masked'].value_counts(normalize=False) 

512 if True in valCounts: 

513 clippedCount = valCounts[True] 

514 else: 

515 clippedCount = 0 

516 print(f'{clippedCount} arc lines where clipped in this iteration of fitting a global dispersion map') 

517 

518 # REMOVE FILTERED ROWS FROM DATA FRAME 

519 mask = (orderPixelTable['residuals_masked'] == True) 

520 orderPixelTable.drop(index=orderPixelTable[ 

521 mask].index, inplace=True) 

522 

523 # a = plt.figure(figsize=(40, 15)) 

524 if arm == "UVB": 

525 fig = plt.figure(figsize=(6, 13.5), constrained_layout=True) 

526 else: 

527 fig = plt.figure(figsize=(6, 11), constrained_layout=True) 

528 gs = fig.add_gridspec(6, 4) 

529 

530 # CREATE THE GRID OF AXES 

531 toprow = fig.add_subplot(gs[0:2, :]) 

532 midrow = fig.add_subplot(gs[2:4, :]) 

533 bottomleft = fig.add_subplot(gs[4:, 0:2]) 

534 bottomright = fig.add_subplot(gs[4:, 2:]) 

535 

536 # ROTATE THE IMAGE FOR BETTER LAYOUT 

537 rotatedImg = np.flipud(np.rot90(self.pinholeFrame.data, 1)) 

538 toprow.imshow(rotatedImg, vmin=10, vmax=50, cmap='gray', alpha=0.5) 

539 toprow.set_title( 

540 "observed arc-line positions (post-clipping)", fontsize=10) 

541 

542 x = orderPixelTable["observed_x"] 

543 # x = np.ones(orderPixelTable.shape[0]) * \ 

544 # self.pinholeFrame.data.shape[1] - orderPixelTable["observed_x"] 

545 toprow.scatter(orderPixelTable["observed_y"], 

546 x, marker='x', c='red', s=4) 

547 

548 # toprow.set_yticklabels([]) 

549 # toprow.set_xticklabels([]) 

550 toprow.set_ylabel("x-axis", fontsize=8) 

551 toprow.set_xlabel("y-axis", fontsize=8) 

552 toprow.tick_params(axis='both', which='major', labelsize=9) 

553 toprow.invert_yaxis() 

554 

555 midrow.imshow(rotatedImg, vmin=10, vmax=50, cmap='gray', alpha=0.5) 

556 midrow.set_title( 

557 "global dispersion solution", fontsize=10) 

558 

559 xfit = orderPixelTable["fit_x"] 

560 # xfit = np.ones(orderPixelTable.shape[0]) * \ 

561 # self.pinholeFrame.data.shape[1] - orderPixelTable["fit_x"] 

562 midrow.scatter(orderPixelTable["fit_y"], 

563 xfit, marker='x', c='blue', s=4) 

564 

565 # midrow.set_yticklabels([]) 

566 # midrow.set_xticklabels([]) 

567 midrow.set_ylabel("x-axis", fontsize=8) 

568 midrow.set_xlabel("y-axis", fontsize=8) 

569 midrow.tick_params(axis='both', which='major', labelsize=9) 

570 midrow.invert_yaxis() 

571 

572 # PLOT THE FINAL RESULTS: 

573 plt.subplots_adjust(top=0.92) 

574 bottomleft.scatter(orderPixelTable["residuals_x"], orderPixelTable[ 

575 "residuals_y"], alpha=0.4) 

576 bottomleft.set_xlabel('x residual') 

577 bottomleft.set_ylabel('y residual') 

578 bottomleft.tick_params(axis='both', which='major', labelsize=9) 

579 

580 hist(orderPixelTable["residuals_xy"], bins='scott', ax=bottomright, histtype='stepfilled', 

581 alpha=0.7, density=True) 

582 bottomright.set_xlabel('xy residual') 

583 bottomright.tick_params(axis='both', which='major', labelsize=9) 

584 subtitle = f"mean res: {mean_res:2.2f} pix, res stdev: {std_res:2.2f}" 

585 fig.suptitle(f"residuals of global dispersion solution fitting - single pinhole\n{subtitle}", fontsize=12) 

586 

587 # GET FILENAME FOR THE RESIDUAL PLOT 

588 filename = filenamer( 

589 log=self.log, 

590 frame=self.pinholeFrame, 

591 settings=self.settings 

592 ) 

593 if self.firstGuessMap: 

594 filename = filename.split("ARC")[0] + "2D_MAP_RESIDUALS.pdf" 

595 else: 

596 filename = filename.split("ARC")[0] + "DISP_MAP_RESIDUALS.pdf" 

597 

598 # plt.show() 

599 home = expanduser("~") 

600 outDir = self.settings["intermediate-data-root"].replace("~", home) 

601 filePath = f"{outDir}/{filename}" 

602 plt.savefig(filePath) 

603 

604 print(f'\nThe dispersion maps fitted against the observed arc-line positions with a mean residual of {mean_res:2.2f} pixels (stdev = {std_res:2.2f} pixles)') 

605 

606 self.log.debug('completed the ``fit_polynomials`` method') 

607 return xcoeff, ycoeff 

608 

609 # use the tab-trigger below for new method 

610 # xt-class-method