Coverage for soxspipe/commonutils/detect_continuum.py : 84%

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*find and fit the continuum in a pinhole flat frame with low-order polynomials. These polynominals are the central loctions of the orders*
6:Author:
7 David Young & Marco Landoni
9:Date Created:
10 September 10, 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 os.path import expanduser
21import numpy as np
22from soxspipe.commonutils.polynomials import chebyshev_xy_polynomial
23from soxspipe.commonutils.filenamer import filenamer
24import matplotlib.pyplot as plt
25from astropy.stats import mad_std
26from astropy.modeling import models, fitting
27from scipy.signal import find_peaks
28from random import random
29from scipy.optimize import curve_fit
30from astropy.stats import sigma_clip, mad_std
31from astropy.visualization import hist
32import collections
33from fundamentals.renderer import list_of_dictionaries
34from soxspipe.commonutils.dispersion_map_to_pixel_arrays import dispersion_map_to_pixel_arrays
35import pandas as pd
38class detect_continuum(object):
39 """
40 *find and fit the continuum in a pinhole flat frame with low-order polynomials. These polynominals are the central loctions of the orders*
42 **Key Arguments:**
43 - ``log`` -- logger
44 - ``pinholeFlat`` -- calibrationed pinhole flat frame (CCDObject)
45 - ``dispersion_map`` -- path to dispersion map csv file containing polynomial fits of the dispersion solution for the frame
46 - ``settings`` -- the settings dictionary
48 **Usage:**
50 To use the ``detect_continuum`` object, use the following:
52 ```python
53 from soxspipe.commonutils import detect_continuum
54 detector = detect_continuum(
55 log=log,
56 pinholeFlat=pinholeFlat,
57 dispersion_map=dispersion_map,
58 settings=settings
59 )
60 order_table_path = detector.get()
61 ```
63 """
65 def __init__(
66 self,
67 log,
68 pinholeFlat,
69 dispersion_map,
70 settings=False
71 ):
72 self.log = log
73 log.debug("instansiating a new 'detect_continuum' object")
74 self.settings = settings
75 self.pinholeFlat = pinholeFlat
76 self.dispersion_map = dispersion_map
78 # KEYWORD LOOKUP OBJECT - LOOKUP KEYWORD FROM DICTIONARY IN RESOURCES
79 # FOLDER
80 self.kw = keyword_lookup(
81 log=self.log,
82 settings=self.settings
83 ).get
84 self.arm = pinholeFlat.header[self.kw("SEQ_ARM")]
86 # DETECTOR PARAMETERS LOOKUP OBJECT
87 self.detectorParams = detector_lookup(
88 log=log,
89 settings=settings
90 ).get(self.arm)
92 # DEG OF THE POLYNOMIALS TO FIT THE ORDER CENTRE LOCATIONS
93 self.polyDeg = self.settings[
94 "soxs-order-centre"]["poly-deg"]
96 return None
98 def get(self):
99 """
100 *return the order centre table filepath*
102 **Return:**
103 - ``order_table_path`` -- file path to the order centre table giving polynomial coeffs to each order fit
104 """
105 self.log.debug('starting the ``get`` method')
107 arm = self.arm
109 # READ THE SPECTRAL FORMAT TABLE TO DETERMINE THE LIMITS OF THE TRACES
110 orderNums, waveLengthMin, waveLengthMax = self.read_spectral_format()
112 # CONVERT WAVELENGTH TO PIXEL POSTIONS AND RETURN ARRAY OF POSITIONS TO
113 # SAMPLE THE TRACES
114 lineList = self.create_pixel_arrays(
115 orderNums,
116 waveLengthMin,
117 waveLengthMax)
118 # SLICE LENGTH TO SAMPLE TRACES IN THE CROSS-DISPERSION DIRECTION
119 self.sliceLength = self.settings[
120 "soxs-order-centre"]["slice-length"]
121 self.peakSigmaLimit = self.settings[
122 "soxs-order-centre"]["peak-sigma-limit"]
124 # FOR EACH ORDER, FOR EACH PIXEL POSITION SAMPLE, FIT A 1D GAUSSIAN IN
125 # CROSS-DISPERSION DIRECTTION. RETURN PEAK POSTIONS
126 lineList = lineList.apply(
127 self.fit_1d_gaussian_to_slice, axis=1)
128 allLines = len(lineList.index)
129 # DROP ROWS WITH NAN VALUES
130 lineList.dropna(axis='index', how='any',
131 subset=['cont_x'], inplace=True)
132 foundLines = len(lineList.index)
133 percent = 100 * foundLines / allLines
134 print(f"{foundLines} out of {allLines} found ({percent:3.0f}%)")
136 from tabulate import tabulate
137 print(tabulate(lineList, headers='keys', tablefmt='psql'))
139 # guassianPixelArrays = {}
140 # for order, xy_array in pixelArrays.items():
141 # found_xy_array = []
142 # for xy in xy_array:
143 # found_xy = self.fit_1d_gaussian_to_slice(
144 # xy, sliceLength, peakSigmaLimit)
145 # if found_xy:
146 # found_xy_array.append(found_xy)
147 # guassianPixelArrays[order] = found_xy_array
149 # percent = 100 * len(found_xy_array) / len(xy_array)
150 # print(f"{len(found_xy_array)} out of {len(xy_array)} found
151 # ({percent:3.0f}%)")
153 # GROUP RESULTS
154 # GET UNIQUE VALUES IN COLUMN
155 uniqueOrders = lineList['Order'].unique()
156 orderGroups = lineList.groupby(['Order'])
158 orderLoctions = {}
159 allResiduals = []
160 allXfit = []
161 allXcoords = []
162 allYcoords = []
163 for o in uniqueOrders:
164 # GROUP RESULTS
165 orderLineList = orderGroups.get_group(o)
166 # ITERATIVELY FIT THE POLYNOMIAL SOLUTIONS TO THE DATA
167 coeff, residuals, xfit, xcoords, ycoords = self.fit_polynomial(
168 orderLineList=orderLineList
169 )
170 allResiduals.extend(residuals)
171 allXfit.extend(xfit)
172 allXcoords.extend(xcoords)
173 allYcoords.extend(ycoords)
174 orderLoctions[o] = coeff
176 # sys.exit(0)
178 # orderLoctions = {}
179 # allResiduals = []
180 # allXfit = []
181 # allXcoords = []
182 # allYcoords = []
183 # for order, pixelArray in guassianPixelArrays.items():
184 # # ITERATIVELY FIT THE POLYNOMIAL SOLUTIONS TO THE DATA
185 # coeff, residuals, xfit, xcoords, ycoords = self.fit_polynomial(
186 # pixelArray=pixelArray
187 # )
188 # allResiduals.extend(residuals)
189 # allXfit.extend(xfit)
190 # allXcoords.extend(xcoords)
191 # allYcoords.extend(ycoords)
192 # orderLoctions[order] = coeff
194 self.plot_results(
195 allResiduals=allResiduals,
196 allXfit=allXfit,
197 allXcoords=allXcoords,
198 allYcoords=allYcoords,
199 orderLoctions=orderLoctions
200 )
202 mean_res = np.mean(np.abs(allResiduals))
203 std_res = np.std(np.abs(allResiduals))
205 print(f'\nThe order centre polynomial fitted against the observed 1D gaussian peak positions with a mean residual of {mean_res:2.2f} pixels (stdev = {std_res:2.2f} pixles)')
207 # WRITE OUT THE FITS TO THE ORDER CENTRE TABLE
208 order_table_path = self.write_order_table_to_file(orderLoctions)
210 self.log.debug('completed the ``get`` method')
211 return order_table_path
213 def read_spectral_format(
214 self):
215 """*read the spectral format table to get some key parameters*
217 **Return:**
218 - ``orderNums`` -- a list of the order numbers
219 - ``waveLengthMin`` -- a list of the maximum wavelengths reached by each order
220 - ``waveLengthMax`` -- a list of the minimum wavelengths reached by each order
221 """
222 self.log.debug('starting the ``read_spectral_format`` method')
224 kw = self.kw
225 pinholeFlat = self.pinholeFlat
226 dp = self.detectorParams
228 # READ THE SPECTRAL FORMAT TABLE FILE
229 home = expanduser("~")
230 calibrationRootPath = self.settings[
231 "calibration-data-root"].replace("~", home)
232 spectralFormatFile = calibrationRootPath + \
233 "/" + dp["spectral format table"]
234 spectralFormat = np.genfromtxt(
235 spectralFormatFile, delimiter=',', names=True)
237 # print(spectralFormat.dtype.names)
239 # EXTRACT REQUIRED PARAMETERS
240 orderNums = spectralFormat["ORDER"]
241 waveLengthMin = spectralFormat["WLMINFUL"]
242 waveLengthMax = spectralFormat["WLMAXFUL"]
244 self.log.debug('completed the ``read_spectral_format`` method')
245 return orderNums, waveLengthMin, waveLengthMax
247 def create_pixel_arrays(
248 self,
249 orderNums,
250 waveLengthMin,
251 waveLengthMax):
252 """*create a pixel array for the approximate centre of each order*
254 **Key Arguments:**
255 - ``orderNums`` -- a list of the order numbers
256 - ``waveLengthMin`` -- a list of the maximum wavelengths reached by each order
257 - ``waveLengthMax`` -- a list of the minimum wavelengths reached by each order
259 **Return:**
260 - ``lineList`` -- a data-frame containing lines and associated pixel locations
261 """
262 self.log.debug('starting the ``create_pixel_arrays`` method')
264 # READ ORDER SAMPLING RESOLUTION FROM SETTINGS
265 sampleCount = self.settings[
266 "soxs-order-centre"]["order-sample-count"]
268 # CREATE THE WAVELENGTH/ORDER ARRAYS TO BE CONVERTED TO PIXELS
269 myDict = {
270 "Order": np.asarray([]),
271 "Wavelength": np.asarray([]),
272 "slit_position": np.asarray([])
273 }
274 for o, wmin, wmax in zip(orderNums, waveLengthMin, waveLengthMax):
276 wlArray = np.arange(
277 wmin, wmax, (wmax - wmin) / sampleCount)
278 myDict["Wavelength"] = np.append(myDict["Wavelength"], wlArray)
279 myDict["Order"] = np.append(
280 myDict["Order"], np.ones(len(wlArray)) * o)
281 myDict["slit_position"] = np.append(
282 myDict["slit_position"], np.zeros(len(wlArray)))
284 lineList = pd.DataFrame(myDict)
285 lineList = dispersion_map_to_pixel_arrays(
286 log=self.log,
287 dispersionMapPath=self.dispersion_map,
288 lineList=lineList
289 )
291 self.log.debug('completed the ``create_pixel_arrays`` method')
292 return lineList
294 def fit_1d_gaussian_to_slice(
295 self,
296 linePixelPostion):
297 """*cut a slice from the pinhole flat along the cross-dispersion direction centred on pixel position, fit 1D gaussian and return the peak pixel position*
299 **Key Arguments:**
300 - ``linePixelPostion`` -- the x,y pixel coordinate from lineList data-frame (series)
302 **Return:**
303 - ``linePixelPostion`` -- now including gaussian fit peak xy position
304 """
305 self.log.debug('starting the ``fit_1d_gaussian_to_slice`` method')
307 # CLIP OUT A SLICE TO INSPECT CENTRED AT POSITION
308 halfSlice = self.sliceLength / 2
310 try:
311 slice = self.pinholeFlat.data[int(linePixelPostion["fit_y"]), max(
312 0, int(linePixelPostion["fit_x"] - halfSlice)):min(2048, int(linePixelPostion["fit_x"] + halfSlice))]
313 except:
314 linePixelPostion["cont_x"] = np.nan
315 linePixelPostion["cont_y"] = np.nan
316 return linePixelPostion
318 # CHECK THE SLICE POINTS IF NEEDED
319 if 1 == 0:
320 x = np.arange(0, len(slice))
321 plt.figure(figsize=(8, 5))
322 plt.plot(x, slice, 'ko')
323 plt.xlabel('Position')
324 plt.ylabel('Flux')
325 plt.show()
327 # EVALUATING THE MEAN AND STD-DEV FOR PEAK FINDING - REMOVES SLICE
328 # CONTAINING JUST NOISE
329 median_r = np.median(slice)
330 std_r = mad_std(slice)
331 peaks, _ = find_peaks(slice, height=median_r +
332 self.peakSigmaLimit * std_r, width=1)
334 # CHECK PEAK HAS BEEN FOUND
335 if peaks is None or len(peaks) <= 0:
336 # CHECK THE SLICE POINTS IF NEEDED
337 if 1 == 0:
338 x = np.arange(0, len(slice))
339 plt.figure(figsize=(8, 5))
340 plt.plot(x, slice, 'ko')
341 plt.xlabel('Position')
342 plt.ylabel('Flux')
343 plt.show()
344 return None
346 # FIT THE DATA USING A 1D GAUSSIAN - USING astropy.modeling
347 # CENTRE THE GAUSSIAN ON THE PEAK
348 g_init = models.Gaussian1D(
349 amplitude=1000., mean=peaks[0], stddev=1.)
350 # print(f"g_init: {g_init}")
351 fit_g = fitting.LevMarLSQFitter()
353 # NOW FIT
354 g = fit_g(g_init, np.arange(0, len(slice)), slice)
355 linePixelPostion["cont_x"] = g.mean + \
356 max(0, int(linePixelPostion["fit_x"] - halfSlice))
357 linePixelPostion["cont_y"] = linePixelPostion["fit_y"]
359 # PRINT A FEW PLOTS IF NEEDED - GUASSIAN FIT OVERLAYED
360 if 1 == 0 and random() < 0.02:
361 x = np.arange(0, len(slice))
362 plt.figure(figsize=(8, 5))
363 plt.plot(x, slice, 'ko')
364 plt.xlabel('Position')
365 plt.ylabel('Flux')
366 guassx = np.arange(0, max(x), 0.05)
367 plt.plot(guassx, g(guassx), label='Gaussian')
368 plt.show()
370 self.log.debug('completed the ``fit_1d_gaussian_to_slice`` method')
371 return linePixelPostion
373 def fit_polynomial(
374 self,
375 orderLineList):
376 """*iteratively fit the dispersion map polynomials to the data, clipping residuals with each iteration*
378 **Key Arguments:**
379 - ``orderLineList`` -- data-frame group containing x,y pixels of measured 1D guassian peak positions
381 **Return:**
382 - ``coeffs`` -- the coefficients of the polynomial fit
383 - ``residuals`` -- the residuals of the fit compared to original guassian peak positions
384 - ``xfit`` -- the polynomial fits x-positions
385 - ``xcoords`` -- the clean guassian peak x-coordinate list (post-clipping)
386 - ``ycoords`` -- the clean guassian peak x-coordinate list (post-clipping)
387 """
388 self.log.debug('starting the ``fit_polynomial`` method')
390 arm = self.arm
392 clippedCount = 1
394 poly = chebyshev_xy_polynomial(
395 log=self.log, deg=self.polyDeg).poly
397 clippingSigma = self.settings[
398 "soxs-order-centre"]["poly-fitting-residual-clipping-sigma"]
400 clippingIterationLimit = self.settings[
401 "soxs-order-centre"]["clipping-iteration-limit"]
403 xcoords = orderLineList["cont_x"].values
404 ycoords = orderLineList["cont_y"].values
406 iteration = 0
407 while clippedCount > 0 and iteration < clippingIterationLimit:
408 iteration += 1
409 # USE LEAST-SQUARED CURVE FIT TO FIT CHEBY POLY
410 coeff = np.ones((self.polyDeg + 1))
411 coeff, pcov_x = curve_fit(
412 poly, xdata=ycoords, ydata=xcoords, p0=coeff)
414 residuals, mean_res, std_res, median_res, xfit = self.calculate_residuals(
415 xcoords=xcoords,
416 ycoords=ycoords,
417 coeff=coeff)
419 # SIGMA-CLIP THE DATA
420 masked_residuals = sigma_clip(
421 residuals, sigma_lower=clippingSigma, sigma_upper=clippingSigma, maxiters=1, cenfunc='median', stdfunc=mad_std)
423 # MASK DATA ARRAYS WITH CLIPPED RESIDUAL MASK
424 startCount = len(xcoords)
425 a = [xcoords, ycoords, residuals]
426 xcoords, ycoords, residuals = [np.ma.compressed(np.ma.masked_array(
427 i, masked_residuals.mask)) for i in a]
428 clippedCount = startCount - len(xcoords)
429 print(f'{clippedCount} pixel positions where clipped in this iteration of fitting an order centre polynomial')
431 self.log.debug('completed the ``fit_polynomials`` method')
432 return coeff, residuals, xfit, xcoords, ycoords
434 def fit_polynomial_bk(
435 self,
436 pixelArray):
437 """*iteratively fit the dispersion map polynomials to the data, clipping residuals with each iteration*
439 **Key Arguments:**
440 - ``pixelArray`` -- the array of x,y pixels of measured 1D guassian peak positions
442 **Return:**
443 - ``coeffs`` -- the coefficients of the polynomial fit
444 - ``residuals`` -- the residuals of the fit compared to original guassian peak positions
445 - ``xfit`` -- the polynomial fits x-positions
446 - ``xcoords`` -- the clean guassian peak x-coordinate list (post-clipping)
447 - ``ycoords`` -- the clean guassian peak x-coordinate list (post-clipping)
448 """
449 self.log.debug('starting the ``fit_polynomial`` method')
451 arm = self.arm
453 clippedCount = 1
455 poly = chebyshev_xy_polynomial(
456 log=self.log, deg=self.polyDeg).poly
458 clippingSigma = self.settings[
459 "soxs-order-centre"]["poly-fitting-residual-clipping-sigma"]
461 clippingIterationLimit = self.settings[
462 "soxs-order-centre"]["clipping-iteration-limit"]
464 xcoords = [v[0] for v in pixelArray]
465 ycoords = [v[1] for v in pixelArray]
467 iteration = 0
468 while clippedCount > 0 and iteration < clippingIterationLimit:
469 iteration += 1
470 # USE LEAST-SQUARED CURVE FIT TO FIT CHEBY POLY
471 coeff = np.ones((self.polyDeg + 1))
472 coeff, pcov_x = curve_fit(
473 poly, xdata=ycoords, ydata=xcoords, p0=coeff)
475 residuals, mean_res, std_res, median_res, xfit = self.calculate_residuals(
476 xcoords=xcoords,
477 ycoords=ycoords,
478 coeff=coeff)
480 # SIGMA-CLIP THE DATA
481 masked_residuals = sigma_clip(
482 residuals, sigma_lower=clippingSigma, sigma_upper=clippingSigma, maxiters=1, cenfunc='median', stdfunc=mad_std)
484 # MASK DATA ARRAYS WITH CLIPPED RESIDUAL MASK
485 startCount = len(xcoords)
486 a = [xcoords, ycoords, residuals]
487 xcoords, ycoords, residuals = [np.ma.compressed(np.ma.masked_array(
488 i, masked_residuals.mask)) for i in a]
489 clippedCount = startCount - len(xcoords)
490 print(f'{clippedCount} pixel positions where clipped in this iteration of fitting an order centre polynomial')
492 self.log.debug('completed the ``fit_polynomials`` method')
493 return coeff, residuals, xfit, xcoords, ycoords
495 def calculate_residuals(
496 self,
497 xcoords,
498 ycoords,
499 coeff):
500 """*calculate residuals of the polynomial fits against the observed line postions*
502 **Key Arguments:**
503 - ``xcoords`` -- the measured x positions of the gaussian peaks
504 - ``ycoords`` -- the measurd y positions of the gaussian peaks
505 - ``coeff`` -- the coefficients of the fitted polynomial
507 **Return:**
508 - ``residuals`` -- x residuals
509 - ``mean`` -- the mean of the residuals
510 - ``std`` -- the stdev of the residuals
511 - ``median`` -- the median of the residuals
512 """
513 self.log.debug('starting the ``calculate_residuals`` method')
515 arm = self.arm
517 poly = chebyshev_xy_polynomial(
518 log=self.log, deg=self.polyDeg).poly
520 # CALCULATE RESIDUALS BETWEEN GAUSSIAN PEAK LINE POSITIONS AND POLY
521 # FITTED POSITIONS
522 xfit = poly(
523 ycoords, *coeff)
524 res_x = np.asarray(xfit) - np.asarray(xcoords)
526 # CALCULATE COMBINED RESIDUALS AND STATS
527 res_mean = np.mean(res_x)
528 res_std = np.std(res_x)
529 res_median = np.median(res_x)
531 self.log.debug('completed the ``calculate_residuals`` method')
532 return res_x, res_mean, res_std, res_median, xfit
534 def plot_results(
535 self,
536 allResiduals,
537 allXfit,
538 allXcoords,
539 allYcoords,
540 orderLoctions):
541 """*generate a plot of the polynomial fits and residuals*
543 **Key Arguments:**
544 - ``allResiduals`` -- list of all residuals
545 - ``allXfit`` -- list of all fitted x-positions
546 - ``allXcoords`` -- cleaned list of all guassian x-pixel positions
547 - ``allYcoords`` -- cleaned list of all guassian y-pixel positions
548 - ``orderLoctions`` -- dictionary of order-location polynomial coeff
550 **Return:**
551 - ``filePath`` -- path to the plot pdf
552 """
553 self.log.debug('starting the ``plot_results`` method')
555 arm = self.arm
557 # a = plt.figure(figsize=(40, 15))
558 if arm == "UVB":
559 fig = plt.figure(figsize=(6, 13.5), constrained_layout=True)
560 else:
561 fig = plt.figure(figsize=(6, 11), constrained_layout=True)
562 gs = fig.add_gridspec(6, 4)
564 # CREATE THE GID OF AXES
565 toprow = fig.add_subplot(gs[0:2, :])
566 midrow = fig.add_subplot(gs[2:4, :])
567 bottomleft = fig.add_subplot(gs[4:, 0:2])
568 bottomright = fig.add_subplot(gs[4:, 2:])
570 # ROTATE THE IMAGE FOR BETTER LAYOUT
571 rotatedImg = np.rot90(self.pinholeFlat.data, 1)
572 toprow.imshow(rotatedImg, vmin=10, vmax=50, cmap='gray', alpha=0.5)
573 toprow.set_title(
574 "1D guassian peak positions (post-clipping)", fontsize=10)
575 x = np.ones(len(allXcoords)) * \
576 self.pinholeFlat.data.shape[1] - allXcoords
577 toprow.scatter(allYcoords, x, marker='x', c='red', s=4)
578 # toprow.set_yticklabels([])
579 # toprow.set_xticklabels([])
580 toprow.set_ylabel("x-axis", fontsize=8)
581 toprow.set_xlabel("y-axis", fontsize=8)
582 toprow.tick_params(axis='both', which='major', labelsize=9)
584 midrow.imshow(rotatedImg, vmin=10, vmax=50, cmap='gray', alpha=0.5)
585 midrow.set_title(
586 "order-location fit solutions", fontsize=10)
587 ylinelist = np.arange(0, self.pinholeFlat.data.shape[0], 3)
588 poly = chebyshev_xy_polynomial(
589 log=self.log, deg=self.polyDeg).poly
590 for o, coeff in orderLoctions.items():
591 xfit = poly(ylinelist, *coeff)
592 xfit = np.ones(len(xfit)) * \
593 self.pinholeFlat.data.shape[1] - xfit
594 xfit, ylinelist = zip(
595 *[(x, y) for x, y in zip(xfit, ylinelist) if x > 0 and x < (self.pinholeFlat.data.shape[1]) - 10])
596 midrow.plot(ylinelist, xfit)
598 # xfit = np.ones(len(xfit)) * \
599 # self.pinholeFrame.data.shape[1] - xfit
600 # midrow.scatter(yfit, xfit, marker='x', c='blue', s=4)
601 # midrow.set_yticklabels([])
602 # midrow.set_xticklabels([])
603 midrow.set_ylabel("x-axis", fontsize=8)
604 midrow.set_xlabel("y-axis", fontsize=8)
605 midrow.tick_params(axis='both', which='major', labelsize=9)
607 # PLOT THE FINAL RESULTS:
608 plt.subplots_adjust(top=0.92)
609 bottomleft.scatter(allXcoords, allResiduals, alpha=0.2, s=1)
610 bottomleft.set_xlabel('x pixel position')
611 bottomleft.set_ylabel('x residual')
612 bottomleft.tick_params(axis='both', which='major', labelsize=9)
614 # PLOT THE FINAL RESULTS:
615 plt.subplots_adjust(top=0.92)
616 bottomright.scatter(allYcoords, allResiduals, alpha=0.2, s=1)
617 bottomright.set_xlabel('y pixel position')
618 bottomright.tick_params(axis='both', which='major', labelsize=9)
619 # bottomright.set_ylabel('x residual')
620 bottomright.set_yticklabels([])
622 mean_res = np.mean(np.abs(allResiduals))
623 std_res = np.std(np.abs(allResiduals))
625 subtitle = f"mean res: {mean_res:2.2f} pix, res stdev: {std_res:2.2f}"
626 fig.suptitle(f"traces of order-centre locations - pinhole flat-frame\n{subtitle}", fontsize=12)
628 # plt.show()
629 filename = filenamer(
630 log=self.log,
631 frame=self.pinholeFlat,
632 settings=self.settings
633 )
634 filename = filename.split("FLAT")[0] + "ORDER_CENTRES_residuals.pdf"
636 home = expanduser("~")
637 outDir = self.settings["intermediate-data-root"].replace("~", home)
638 filePath = f"{outDir}/{filename}"
639 plt.savefig(filePath)
641 self.log.debug('completed the ``plot_results`` method')
642 return filePath
644 def write_order_table_to_file(
645 self,
646 orderLoctions):
647 """*write out the fitted polynomial solution coefficients to file*
649 **Key Arguments:**
650 - ``orderLoctions`` -- dictionary of the order coefficients
652 **Return:**
653 - ``order_table_path`` -- path to the order table file
654 """
655 self.log.debug('starting the ``write_order_table_to_file`` method')
657 arm = self.arm
659 # SORT COEFFICIENT OUTPUT TO WRITE TO FILE
660 listOfDictionaries = []
661 for k, v in orderLoctions.items():
662 orderDict = collections.OrderedDict(sorted({}.items()))
663 orderDict["order"] = k
664 orderDict["degy"] = self.polyDeg
665 n_coeff = 0
666 for i in range(0, self.polyDeg + 1):
667 orderDict[f'CENT_c{i}'] = v[n_coeff]
668 n_coeff += 1
669 listOfDictionaries.append(orderDict)
671 # DETERMINE WHERE TO WRITE THE FILE
672 home = expanduser("~")
673 outDir = self.settings["intermediate-data-root"].replace("~", home)
675 filename = filenamer(
676 log=self.log,
677 frame=self.pinholeFlat,
678 settings=self.settings
679 )
680 filename = filename.split("FLAT")[0] + "ORDER_CENTRES.csv"
682 order_table_path = f"{outDir}/{filename}"
683 dataSet = list_of_dictionaries(
684 log=self.log,
685 listOfDictionaries=listOfDictionaries
686 )
687 csvData = dataSet.csv(filepath=order_table_path)
689 self.log.debug('completed the ``write_order_table_to_file`` method')
690 return order_table_path
692 # use the tab-trigger below for new method
693 # xt-class-method