Skip to content
Snippets Groups Projects

Plotten eines Intervalls um eine Spektrallinie herum (aus FITS-Datei)

  • Clone with SSH
  • Clone with HTTPS
  • Embed
  • Share
    The snippet can be accessed without any authentication.
    Authored by Marco Ammon

    Soll man für das Spektroskopie-Protokoll etwa eine bestimmte Spektrallinie plotten, kann man dieses Tool verwenden.

    Benutzung:

    python3 interval.py FITS-DATEI SPEKTRALLINIE BEREICH AUSGABEDATEI

    • FITS-DATEI gibt das zu plottende Spektrum an, etwa stern3.fit

    • SPEKTRALLINIE gibt die darzustellende Spektrallinie in Angstrom an, etwa 4685.68

    • BEREICH gibt die den darzustellenden Bereich in Angstrom an (also von SPEKTRALLINIE - BEREICH bis SPEKTRALLINIE + BEREICH), etwa 4 um einen 8 Angstrom breiten Plot zu erhalten

    • AUSGABEDATEI gibt den Pfad der Zieldatei an, etwa stern3-he2.pdf

    Ergebnis:

    He 2-Linie bei Stern3

    Edited
    interval.py 1.12 KiB
    from astropy.io import fits
    import numpy as np
    import matplotlib.pyplot as plt
    import sys
    
    def main():
        with fits.open(sys.argv[1]) as hdul:
            hdul.info()
            hdu = hdul[0]
            hdu.verify('fix')
            data = hdu.data
            header = hdu.header
            plot_intensity(header, data)
            plt.tight_layout()
            plt.savefig(sys.argv[4])
    
    
    def plot_intensity(header, data):
        wavelength_of_interest = float(sys.argv[2])
        interval_of_interest = float(sys.argv[3])
        start_wavelength = header['CRVAL1']
        increment_wavelength = header['CDELT1']
        pixel_number = len(data)
        wavelengths = start_wavelength + np.arange(0, pixel_number) * increment_wavelength
        condition = np.logical_and(wavelengths >= (wavelength_of_interest - interval_of_interest), wavelengths <= (wavelength_of_interest + interval_of_interest))
        wavelengths = np.extract(condition, wavelengths)
        intensities = np.extract(condition, data)
        plt.gcf().set_size_inches(7,4.5)
        plt.xlabel('Wellenlänge $[\mathrm{\AA}]$')
        plt.ylabel('Intensität')
        plt.plot(wavelengths, intensities)
    
    if __name__ == "__main__":
        main()
    
    0% Loading or .
    You are about to add 0 people to the discussion. Proceed with caution.
    Finish editing this message first!
    Please register or to comment