GrisFitsFile.py 19.2 KB
Newer Older
Carl Schaffer's avatar
Carl Schaffer committed
1
2
3
4
5
6
"""
Created on 2018-04-09
@author Carl Schaffer
@mail   carl.schaffer@leibniz-kis.de
"""

Carl Schaffer's avatar
Carl Schaffer committed
7
import datetime  # date handling
8
import logging
9
10
import os  # path and file manipulations
import sys  # error handling
11
from glob import glob
12
from os.path import basename, exists, join, dirname
13
from pathlib import Path
14

Carl Schaffer's avatar
linting    
Carl Schaffer committed
15
import astropy.io.fits as fitsio  # fits reading and writing
16
import astropy.units as u
17
import matplotlib.pyplot as plt
Carl Schaffer's avatar
Carl Schaffer committed
18
import numpy as np
19
20
21
import pytz
from astroplan import Observer
from astropy.coordinates import SkyCoord, Angle, get_sun
22
from astropy.wcs import WCS
Carl Schaffer's avatar
Carl Schaffer committed
23

Carl Schaffer's avatar
Carl Schaffer committed
24
from kis_tools.generic.fits import FitsFile
25
from kis_tools.gris.coord_correction_ml import get_coords_ml
Carl Schaffer's avatar
Carl Schaffer committed
26
from kis_tools.gris.util import telescope_to_hpc, expand_coords, calculate_slit_translation
27
from kis_tools.util.calculations import estimate_telescope_drift
28
from kis_tools.util.constants import gregor_coords, gris_stepsize_along_slit
Carl Schaffer's avatar
Carl Schaffer committed
29
from kis_tools.util.util import date_from_fn, gris_obs_mode
Carl Schaffer's avatar
Carl Schaffer committed
30

31
32
logger = logging.getLogger(__name__)

Carl Schaffer's avatar
Carl Schaffer committed
33

34
class GrisFitsFile(FitsFile):
Carl Schaffer's avatar
Carl Schaffer committed
35
36
37
38
39
    """Class for handling fits files for GRIS db
    When the constructor is called, the filename is stored and basic
    metadata parameters such as date and runnumber are extracted from
    the filename. Once the class function GrisFitsFile.parse() is
    called, metadata from the fits file is available as a dictionary
Carl Schaffer's avatar
typo    
Carl Schaffer committed
40
    stored at GrisFitsFile.properties. Any errors occurring during
Carl Schaffer's avatar
Carl Schaffer committed
41
42
    opening and parsing of the file are stored as strings in
    GrisFitsFile.errors
Carl Schaffer's avatar
Carl Schaffer committed
43

44
45
46
47
48
    Args:
      filename: filename of target fits file

    Returns:

Carl Schaffer's avatar
Carl Schaffer committed
49
    """
Carl Schaffer's avatar
Carl Schaffer committed
50

Carl Schaffer's avatar
Carl Schaffer committed
51
52
    def __init__(self, filename):
        """Constructor"""
53
54
        if isinstance(filename, Path):
            filename = str(filename.resolve())
55
        super().__init__(filename)
Carl Schaffer's avatar
Carl Schaffer committed
56
        assert exists(filename)
Carl Schaffer's avatar
Carl Schaffer committed
57
58
59
60
        self.properties = {}
        self.errors = []

        self.filename = filename
61
        self.name = basename(filename)
62
        self.path = self.filename
63

64
        self.dirname = dirname(self.path)
Carl Schaffer's avatar
Carl Schaffer committed
65
66
67
68
        self.parse_filename()
        self.is_parsed = False

    def parse_filename(self):
69
        """Extract metadata from filename"""
Carl Schaffer's avatar
Carl Schaffer committed
70
        # Extract time from Filename
71
        name = os.path.splitext(basename(self.filename))[0].split("_")
Carl Schaffer's avatar
linting    
Carl Schaffer committed
72
73
74
        year = int(name[1][:4])
        month = int(name[1][4:6])
        day = int(name[1][6:])
Carl Schaffer's avatar
Carl Schaffer committed
75
76

        # Get time
Carl Schaffer's avatar
linting    
Carl Schaffer committed
77
78
79
        # hh = int(name[2][:2])
        # mm = int(name[2][2:4])
        # outstring = int(name[2][4:6])
Carl Schaffer's avatar
Carl Schaffer committed
80

Carl Schaffer's avatar
linting    
Carl Schaffer committed
81
        # self.date = datetime.datetime(year,month,day,hh,mm,outstring)
Carl Schaffer's avatar
Carl Schaffer committed
82
        self.date = datetime.datetime(year, month, day, 0, 0, 0)
Carl Schaffer's avatar
Carl Schaffer committed
83

84
        # Extract Run-, Map-, and Stepnumber, as well as mode string
85
        self.runnumber = int(name[4])
86
        self.mapnumber = int(name[5])
87
        self.slitposnumber = int(name[6])
88
89
        self.modestring = name[3]

90
    @property
Carl Schaffer's avatar
Carl Schaffer committed
91
    def verbose_runname(self):
92
93
94
95
96
97
        """
        Get verbose run name such as 24apr14.002

        Returns:
            outstring: verbose run name
        """
98
99
100
101
        outstring = self.date.strftime("%d%b%y").lower()
        outstring += f".{int(self.runnumber):03d}"
        return outstring

102
103
104
105
106
107
108
109
110
111
112
113
114
    @property
    def telescope_centering(self):
        """
        Determine state of telescope centering at the time of observation, queries GDBS
        for centering data

        Returns:
            outstring: x_offset, y_offset, x_uncertainty, y_uncertainty all in arcseconds
        """
        # Query centering from GDBS
        try:
            from kis_tools.gdbs.query_gdbs import query_gdbs
            gdbs_usable = True
115
116
            x_telcenter, x_valid, x_age = query_gdbs('xsuncenter', self.obs_time)
            y_telcenter, y_valid, y_age = query_gdbs('ysuncenter', self.obs_time)
117
118
119
120
121
122
123
124
125
126
            # Determine uncertainties for centering values
            centering_age = max(x_age, y_age)
            drift = estimate_telescope_drift(centering_age.seconds / 3600)
            x_std, y_std = drift, drift

            # if the centering is older than 24 hours, discard
            if centering_age > datetime.timedelta(hours=24):
                gdbs_usable = False
        except IOError:
            gdbs_usable = False
127
128
        except ValueError:
            gdbs_usable = False
129

130
        # Catch instances where no numbers but error strings are written to GDBS
Carl Schaffer's avatar
bugfix    
Carl Schaffer committed
131
132
133
134
135
        if gdbs_usable:
            try:
                x_telcenter, y_telcenter = float(x_telcenter), float(y_telcenter)
            except ValueError:
                gdbs_usable = False
136

137
138
139
140
141
142
143
        if not gdbs_usable:
            # Determined from offset between uncentered and cross_correlated data
            # STD-DEV of distances in X and Y again see gris coord study for details
            uncertainty_no_center_x = 117.1
            uncertainty_no_center_y = 182.4
            x_telcenter, y_telcenter, x_std, y_std = 0, 0, uncertainty_no_center_x, uncertainty_no_center_y

144
        return x_telcenter, y_telcenter, x_std, y_std
145

Carl Schaffer's avatar
Carl Schaffer committed
146
    @property
147
    def _coords_from_simple_header(self):
148
149
150
151
152
153
154
155
156
157
        """coordinates for slit, the Keywords 'SLITPOSX' and 'SLITPOSY' are assumed to be the center of the slit
         while the keyword 'SLITORIE' describes the rotation of the slit w.r.t. the helioprojective axes. The
         algorithm further assumes that each pixel represents a square area with side length of 'STEPSIZE'.
         These assumptions are used to calculate the coordinates for the center of each pixel within the data array.

        Returns:
            (X,Y):ndarray(float) x and y coordinates for each pixel of the slit


        WARNING! The code below is not checked for correctness! We do not know what direction SLITORIE references
Carl Schaffer's avatar
Carl Schaffer committed
158
159
        off and coordinates provided by GREGOR are very error-prone. If you need precision better than ~50" this
        function is not sufficient.
160
161
        """

Carl Schaffer's avatar
Carl Schaffer committed
162
        h = self.header
Carl Schaffer's avatar
Carl Schaffer committed
163
        slit_length_pixels = self.NAXIS2 if self.old_header else self.NAXIS1
Carl Schaffer's avatar
Carl Schaffer committed
164

165
        # Slit center in telescope coordinates
Carl Schaffer's avatar
Carl Schaffer committed
166
        x_slitcenter = self.SLITPOSX if self.header['SLITPOSX'] != "" else 0
Carl Schaffer's avatar
Carl Schaffer committed
167
        y_slitcenter = self.SLITPOSY if self.header['SLITPOSY'] != "" else 0
Carl Schaffer's avatar
Carl Schaffer committed
168
169

        coord_ref = slit_length_pixels // 2
Carl Schaffer's avatar
Carl Schaffer committed
170
171
172
173
        if "SLITCNTR" in self.header:
            if self.SLITCNTR.strip() != 'scancntr':
                raise ValueError(
                    f'Cannot determine reference pixel for coordinates! SLITCNTR is {self.SLITCNTR} instead of scancntr')
Carl Schaffer's avatar
Carl Schaffer committed
174

175
176
177
178
179
        # Transform to HPC (shift by telescope centering, rotate by p0 angle)
        delta_x, delta_y, std_delta_x, std_delta_y = self.telescope_centering

        # If uncertainties are too large, use ML model to correct the coordinates
        from kis_tools.gris.coord_correction_ml import coord_model_stdx, coord_model_stdy
180

181
        if std_delta_x > coord_model_stdx or std_delta_y > coord_model_stdy:
182
183
184
185
186
187
            try:
                x_hpc, y_hpc, std_delta_x, std_delta_y = get_coords_ml(self)
            except ValueError as e:
                logger.warning(f"Error in finding ML coords for {self.__repr__()}:{e}")
                # Point to slitpos with 1 solar disk uncertainty if coord retrieval didn't work
                x_hpc, y_hpc, std_delta_x, std_delta_y = x_slitcenter, y_slitcenter, 1000, 1000
188
189
        else:
            x_hpc, y_hpc = telescope_to_hpc(x_slitcenter, y_slitcenter, self.p0_angle, delta_x, delta_y)
Carl Schaffer's avatar
Carl Schaffer committed
190

191
192
193
194
195
        try:
            angle = self.slit_orientation
        except ValueError as e:
            logger.warning(f"Could not determine Slit orientation! \n {e}")
            angle = 0
Carl Schaffer's avatar
Carl Schaffer committed
196

197
198
        X, Y = expand_coords(x_hpc, y_hpc,
                             slit_orientation=angle,
Carl Schaffer's avatar
Carl Schaffer committed
199
                             reference_pixel=coord_ref,
200
201
202
                             stepsize_slit=gris_stepsize_along_slit,
                             slit_length_pixels=slit_length_pixels)

203
204
205
206
207
        coord_shift = calculate_slit_translation(slit_number=self.slitposnumber, stepsize_perp=self.STEPSIZE,
                                                 slit_orientation=angle)

        X = X + coord_shift[0]
        Y = Y + coord_shift[1]
208
209
210
211
212

        coord_array = np.array([X, Y]) * u.arcsec
        # reorder axes for consistency with other coordinate retrieval methods
        coord_array = np.swapaxes(coord_array, 0, 1)

213
214
        # track coordinate uncertainty within class
        self.coord_uncertainty = (std_delta_x, std_delta_y)
215
        return coord_array
216

217
218
    @property
    def coords(self):
Carl Schaffer's avatar
Carl Schaffer committed
219
220
221
222
        """
        Returns:array(float) shape (2,len(slit)) helioprojective coordinates in arcsecnd along slit
        """

223
224
225
226
227
228
        try:
            coords = self._coords_from_wcs
        except AssertionError:
            coords = self._coords_from_simple_header
        return coords

229
230
231
232
    @property
    def old_header(self):
        """Check if the header of this file has already been processed"""
        old_header = True
Carl Schaffer's avatar
Carl Schaffer committed
233
        if any([k in self.header for k in ["SOLARNET"]]):
234
235
            old_header = False
        return old_header
Carl Schaffer's avatar
Carl Schaffer committed
236

237
238
239
240
241
    @property
    def zenit(self):
        """
        Get Zenit angle (90° - elevation)
        Returns: zenit float in degrees
242
        Raises Value error if elevation is not contained in the header
243
244

        """
245
246
247
248
        elev = self._get_value(['ELEV_ANG', 'ELEVATIO'])
        if elev == '':
            raise ValueError(f'No elevation tracked for {self.__repr__()}')
        zenit = 90 - elev
249
250
        return zenit

Carl Schaffer's avatar
Carl Schaffer committed
251
252
    @property
    def slit_orientation(self):
Carl Schaffer's avatar
Carl Schaffer committed
253
254
        """Get slitorintation:
         clockwise angle from the solar equator, scan direction is 90° further in clockwise direction.
255

256
        Attempts to correct for derotator and image rotation, usually correct within a couple 10s of degrees.
257
        Returns:
258
            slitorie: float (degrees)
259
        """
260

Carl Schaffer's avatar
Carl Schaffer committed
261
262
263
        # Pull values from header

        # p0 angle, rotates over the course of the year
264
        p0_angle = self.p0_angle
Carl Schaffer's avatar
Carl Schaffer committed
265
266

        # Azimut and zenit at time of observation start
267
268
        azimut = self._get_value(['AZIMUT'])
        zenit = self.zenit
Carl Schaffer's avatar
Carl Schaffer committed
269
270

        # Slitorientation from header
271
        slitorie = self.SLITORIE
272

Carl Schaffer's avatar
Carl Schaffer committed
273
        # Calculate declination of disk center from time
274
275
276
277
278
        gregor_observer = Observer(location=gregor_coords)
        date_tz = pytz.utc.localize(self.obs_time)
        time = gregor_observer.datetime_to_astropy_time(date_tz)
        sun = get_sun(time)
        dec = sun.dec
Carl Schaffer's avatar
Carl Schaffer committed
279
280
281
282
283
284
285
286

        # Calculate parallactic angle as described in the GREGOR derotator manual
        # http://www.leibniz-kis.de/fileadmin/user_upload/observatorien/gre_docs/man/GRE-KIS-MAN-0008_v001_derotator.pdf
        parallactic_angle = np.degrees(
            np.arcsin(np.sin(np.radians(azimut)) * np.cos(gregor_coords.lat.to(u.rad)) / np.cos(dec.to(u.rad))))
        self.parallactic_angle = parallactic_angle

        const_term = -130  # determined by trial and error to match results from cross correlation with HMI
287
288

        alpha = -1 * parallactic_angle.to(u.degree).value + azimut - zenit - p0_angle + const_term
289
        if self.derotator_in:
290
            rotangle = self.ROTANGLE
291
292
293
294
            alpha = -2 * rotangle + (const_term / 2)
        angle = slitorie - alpha

        return angle
Carl Schaffer's avatar
Carl Schaffer committed
295
296

    @property
297
298
    def number_of_maps(self):
        """
299
300

        Returns:
301
            nmaps: int number of maps in associated observation
302
        """
Carl Schaffer's avatar
Carl Schaffer committed
303
        return self._get_value(["NMAPS", "SERIES"])
Carl Schaffer's avatar
Carl Schaffer committed
304
305
306

    @property
    def map_index(self):
307
308
309
310
311
        """

        Returns:
            i_map: int map index within observation
        """
Carl Schaffer's avatar
Carl Schaffer committed
312
        return self._get_value(["IMAP", "ISERIE"])
Carl Schaffer's avatar
Carl Schaffer committed
313
314
315

    @property
    def n_steps(self):
316
317
318
319
320
321
322
323
324
325
326
327
328
        """
        Determine the number of slit scanning steps in the associated map.

        Returns:
            nsteps: number of steps in the associated map

        """
        # get initial guess from header
        nsteps = self._get_value(["STEPS", "NSTEPS"])

        # check if other files matching the same run are present. If so, check the running step number
        # for the same map to determine the actual number of steps. It is assumed, that all files belonging to the
        # same map are stored in the same directory
Carl Schaffer's avatar
Carl Schaffer committed
329
        pattern = f"gris_*_*_*_{self.runnumber:03d}_{self.mapnumber:03d}_*.fits"
330
331
        matching_files = glob(join(self.dirname, pattern))
        max_steps = int(
Carl Schaffer's avatar
Carl Schaffer committed
332
            max([basename(m).split(".")[0].split("_")[-1] for m in matching_files])
333
334
335
336
337
338
        )

        if max_steps < nsteps:
            nsteps = max_steps

        return nsteps
Carl Schaffer's avatar
Carl Schaffer committed
339

340
341
    @property
    def obs_time(self):
342
343
344
345
346
        """

        Returns:
            obs_time: observation time given in filename
        """
Carl Schaffer's avatar
Carl Schaffer committed
347
        return date_from_fn(self.filename)
348

Carl Schaffer's avatar
Carl Schaffer committed
349
350
    @property
    def step_index(self):
351
352
353
354
355
        """

        Returns:
            step_index: index of scanning step within current map
        """
Carl Schaffer's avatar
Carl Schaffer committed
356
        return self._get_value(["ISTEP"])
Carl Schaffer's avatar
Carl Schaffer committed
357

358
359
360
361
362
363
364
    @property
    def p0_angle(self):
        """

        Returns:
            p0_angle: angle between telescope and solar north
        """
365
366
367
368
        angle = self._get_value(["P0ANGLE", 'SOLAR_P0'])
        if angle == '':
            raise ValueError(f'Header of {self} contains illegal p0_angle:"{angle}"')
        return angle
369

370
371
    @property
    def step_angle(self):
372
373
374
375
376
377
        """

        Returns:
            step_angle: scanning angle
        """

Carl Schaffer's avatar
Carl Schaffer committed
378
        return self._get_value(["STEPANGL"])
379
380
381
382

    @property
    def wavelength_step_size(self):
        """Get step length along wavelength axis. We use the values from the continuum correction fitting"""
Carl Schaffer's avatar
Carl Schaffer committed
383
        return self._kw_mean(["FF1WLDSP", "FF2WLDSP"])
384
385
386
387
388

    @property
    def wavelength_region(self):
        """Get wavelength region covered by file. Returns tuple of (upper, lower)"""
        step_size = self.wavelength_step_size
389

390
        n_steps_wl = self.data.shape[-1] if self.old_header else self.NAXIS3
391
392

        # average for lower border:
Carl Schaffer's avatar
Carl Schaffer committed
393
394
        lower = self._kw_mean(["FF1WLOFF", "FF2WLOFF"])
        # add wl steps to calculate upper border
395
396
397
398
        upper = lower + step_size * n_steps_wl

        return lower, upper

399
400
401
402
403
404
405
406
407
    @property
    def derotator_in(self):
        """
        Determine whether the derotator is inserted

        Returns:
            bool: derotator in

        """
408
409
410
411
412
413
        try:
            rotval = self._get_value(["ROTCODE"])
        except KeyError:
            # if no derotator keywords are scontained, assume it doesn't exist yet
            return False

414
415
416
417
418
        if rotval == 0:
            return True
        else:
            return False

Carl Schaffer's avatar
Carl Schaffer committed
419
420
421
422
423
424
425
426
427
428
429
430
    @property
    def states(self):
        """
        Determine number of polarimetric states

        Returns:
            n_states int: number of states
        """

        # Number of polarimetric states is stored in the last axis. Determine
        # number of axes from NAXIS keyword and query value for last axis.
        naxes = self.query("NAXIS")
431
        keyword = f"NAXIS{naxes}"
Carl Schaffer's avatar
Carl Schaffer committed
432
433
434
435
        states = self.query(keyword)

        return states

436
437
    @property
    def obs_mode(self):
438
        """Get observation mode"""
439
        self.parse()
440
441
442
443
        # extract obs mode keywords
        states = self.states
        n_maps = self.number_of_maps
        n_steps = self.n_steps
444

445
        mode = gris_obs_mode(states, n_maps, n_steps)
446
        return mode
Carl Schaffer's avatar
Carl Schaffer committed
447

Carl Schaffer's avatar
Carl Schaffer committed
448
449
450
    def parse(self):
        """extract data and metadata from file"""
        if self.is_parsed:
Carl Schaffer's avatar
Carl Schaffer committed
451
            return
452
453
454
455
456
457
458
459
460

        def is_valid(value):
            valid = True
            if isinstance(value, fitsio.card.Undefined):
                valid = False
            elif value == "-NAN":
                valid = False

            return valid
Carl Schaffer's avatar
Carl Schaffer committed
461
462

        # Read file
Carl Schaffer's avatar
Carl Schaffer committed
463
        try:
Carl Schaffer's avatar
Carl Schaffer committed
464
            header = self.header
465

Carl Schaffer's avatar
Carl Schaffer committed
466
            # Add header to self.properties, cast strange data types where
467
468
            # necessary, fill invalid values with empty strings

Carl Schaffer's avatar
Carl Schaffer committed
469
            for key in header.keys():
470
471
472
                value = header[key]
                # check validity:
                value = value if is_valid(value) else ""
Carl Schaffer's avatar
Carl Schaffer committed
473
                if key in ["COMMENT", "HISTORY"]:
474
                    self.properties[key] = str(value)
Carl Schaffer's avatar
Carl Schaffer committed
475
                else:
476
477
                    self.properties[key] = value

Carl Schaffer's avatar
Carl Schaffer committed
478
            self.is_parsed = True
Carl Schaffer's avatar
Carl Schaffer committed
479
        except:
Carl Schaffer's avatar
linting    
Carl Schaffer committed
480
            exception = sys.exc_info()[0]
Carl Schaffer's avatar
Carl Schaffer committed
481
            self.errors.append("File %s" % self.filename + str(exception))
Carl Schaffer's avatar
Carl Schaffer committed
482
483

    def make_bson(self):
Carl Schaffer's avatar
Carl Schaffer committed
484
485
        """construct BSON document for mongodb storage"""
        return self.properties
Carl Schaffer's avatar
Carl Schaffer committed
486

487
    @property
488
    def _coords_from_wcs(self):
489
490
491
492
        """
        Retrieve spatial coordinates for each pixel in for slices projected into the spatial dimensions.
        Returns a data cube shaped s the spatial dimensions
        Returns: arcsec, cube of (NAXIS1, NAXIS2, 2) Values are Helioprojective X and Y in arcseconds
493
494
495

        Side-effect: Sets self.coord_uncertainty

496
497
        """

Carl Schaffer's avatar
Carl Schaffer committed
498
499
500
        assert (
            not self.old_header
        ), f"Can't get coordinates from {self.path} ! Run header conversion first!"
501

502
503
504
        wcs = WCS(self.header)
        pixel_shape = wcs.pixel_shape
        n_y = pixel_shape[1]
505
        slit_list = np.zeros((n_y, self["NAXIS"]))
506
507
        slit_list[:, 1] = np.arange(n_y)
        slit_list[:, 0] = 0
508
        converted = np.array(wcs.pixel_to_world_values(slit_list))
509
        degrees = Angle(u.degree * converted[:, :2]).wrap_at(180 * u.degree)
510
        arcsec = degrees.to(u.arcsec)
511

512
513
514
        # Track coord uncertainty within instance
        self.coord_uncertainty = (self.CSYER1, self.CSYER2)

515
516
        return arcsec

517
    def plot_slit(self):
518
        m = self.full_disk_map
519

520
        # retrieve different data sources
521
        coordinates = self.coords
522
523
524
525

        X, Y = coordinates[:, 0], coordinates[:, 1]

        # determine fov limits for plotting
526

527
        plotting_padding = 300 * u.arcsec
528
529
530
531
        x_min = min(X) - plotting_padding
        x_max = max(X) + plotting_padding
        y_min = min(Y) - plotting_padding
        y_max = max(Y) + plotting_padding
532
533

        # Create submap of FOV and plot
Carl Schaffer's avatar
Carl Schaffer committed
534
535
536
537
        submap = m.submap(
            SkyCoord(x_min, y_min, frame=m.coordinate_frame),
            SkyCoord(x_max, y_max, frame=m.coordinate_frame),
        )
538
539
540
541
542
543
        fig = plt.figure()  # noqa F841
        ax = plt.subplot(projection=submap)
        im = submap.plot()  # noqa:F841

        # create SkyCoord objects and overplot
        slit_wcs = SkyCoord(X, Y, frame=submap.coordinate_frame)
Carl Schaffer's avatar
Carl Schaffer committed
544
        ax.plot_coord(slit_wcs, "g-", label="WCS Coordinates")
545
546
547

        # Create Legend and annotations
        title = "Gris Slit {0}_{1:03d}_{2:03d}_{3:03d}"
Carl Schaffer's avatar
Carl Schaffer committed
548
549
550
551
552
553
        title = title.format(
            self.obs_time.strftime("%Y%m%d"),
            self.runnumber,
            self.mapnumber,
            self.step_index,
        )
554
555
        plt.title(title)

Carl Schaffer's avatar
Carl Schaffer committed
556
557
558
559
560
561
562
563
564
        annotation = (
            f"Spectrum taken at \n{self.obs_time.strftime('%Y-%m-%d %H:%M:%S')}"
        )
        annotation += (
            f"\nHMI Reference taken at \n{submap.date.strftime('%Y-%m-%d %H:%M:%S')}"
        )
        plt.text(
            0.05, 0.8, annotation, horizontalalignment="left", transform=ax.transAxes
        )
565
566
        plt.tight_layout(pad=3.5)
        _ = plt.legend()
567
568
569

        return fig, ax

Carl Schaffer's avatar
Carl Schaffer committed
570
    def __str__(self):
Carl Schaffer's avatar
Carl Schaffer committed
571
572
573
        outstring = (
            "GrisFitsFile %s from %s:\n\tRun: %s\n\tMap: %s\n\tSlitposition: %s\n"
        )
Carl Schaffer's avatar
linting    
Carl Schaffer committed
574
        outstring = outstring % (
Carl Schaffer's avatar
Carl Schaffer committed
575
            os.path.basename(self.filename),
Carl Schaffer's avatar
Carl Schaffer committed
576
577
            self.date,
            self.runnumber,
578
            self.map_index,
Carl Schaffer's avatar
Carl Schaffer committed
579
            self.slitposnumber,
Carl Schaffer's avatar
Carl Schaffer committed
580
        )
Carl Schaffer's avatar
linting    
Carl Schaffer committed
581
        return outstring
582
583


Carl Schaffer's avatar
Carl Schaffer committed
584
if __name__ == "__main__":
585
586
587
    gf = GrisFitsFile(sys.argv[1])
    gf.parse()
    print(gf)