attoworld.spectrum

32def wavelength_to_frequency(
33    wavelengths_nm: np.ndarray,
34    spectrum: np.ndarray,
35    frequencies: Optional[np.ndarray] = None,
36):
37    """Convert a wavelength spectrum in W/nm into a frequency spectrum in W/THz
38
39    Args:
40        wavelength_nm (np.ndarray): the wavelengths included in the data, in nanometers
41        spectrum (np.ndarray): the spectrum corresponding to the input wavelength scale
42        frequencies: (optional) frequency vector for the output. If not specified, a vector will be calculated such that resolution and range are preserved.
43    Returns:
44        f: frequencies (Hz)
45        scaled_spectrum: the frequency-domain spectrum
46    """
47    # Contributed by Nick Karpowicz
48    input_frequencies = 1e9 * constants.speed_of_light / wavelengths_nm
49
50    if frequencies is None:
51        frequency_step = np.min(np.abs(np.diff(input_frequencies)))
52        min_frequency = np.min(input_frequencies)
53        max_frequency = np.max(input_frequencies)
54        frequency_count = int(np.ceil(max_frequency - min_frequency) / frequency_step)
55        frequencies = min_frequency + frequency_step * np.array(
56            range(0, frequency_count), dtype=float
57        )
58
59    # apply Jakobian scaling and use W/THz
60    scaled_spectrum = (
61        1e12 * constants.speed_of_light * spectrum / (input_frequencies**2)
62    )
63
64    scaled_spectrum = interpolate(
65        frequencies, input_frequencies, scaled_spectrum, inputs_are_sorted=False
66    )
67
68    return frequencies, scaled_spectrum

 8def frequency_to_wavelength(
 9    frequencies: np.ndarray,
10    spectrum: np.ndarray,
11    wavelengths: Optional[np.ndarray] = None,
12):
13    """Convert a frequency spectrum in W/Hz into a wavelength spectrum in W/m.
14    SI units.
15    Args:
16        frequencies (np.ndarray): the frequencies included in the data, in Hz
17        spectrum (np.ndarray): the spectrum corresponding to the input frequency scale
18        wavelengths: (optional) wavelength vector for the output. If not specified, a the output data will be on the same grid, but scaled
19    Returns:
20        wavelengths: wavelengths in m
21        scaled_spectrum: the wavelength-domain spectrum
22    """
23    if wavelengths is None:
24        wavelengths = constants.speed_of_light / frequencies[frequencies > 0.0]
25        spectrum = spectrum[frequencies > 0.0]
26    else:
27        wavelengths = wavelengths[wavelengths > 0.0]
28
29    return wavelengths, spectrum / (wavelengths**2)

71def transform_limited_pulse_from_spectrometer(
72    wavelengths_nm: np.ndarray, spectrum: np.ndarray, gate_level: Optional[float] = None
73):
74    """Calculates the transform-limited pulse corresponding to a spectrum
75
76    Args:
77        wavelength_nm: the wavelengths included in the data, in nanometers
78        spectrum: the spectrum corresponding to the input wavelength scale
79        gate_level: (optional) level, relative to the maximum at which to apply a gate to the spectrum. For example, with gate_level=0.01, values less than 1% of the maximum signal will be set to zero
80
81    Returns:
82        t: time vector (s)
83        pulse: the pulse intensity vs. time
84    """
85    f, spec = wavelength_to_frequency(wavelengths_nm, spectrum)
86    df = f[1] - f[0]
87    t = np.fft.fftshift(np.fft.fftfreq(f.shape[0], d=df))
88    gated_spectrum = np.array(spec)
89    if gate_level is not None:
90        gated_spectrum[spec < (gate_level * np.max(spec))] = 0.0
91
92    pulse = np.fft.fftshift(np.abs(np.fft.ifft(np.sqrt(gated_spectrum)))) ** 2
93
94    return t, pulse

10def load_calibration_data(calibration_data_filepath):
11    """
12    Load calibration data in a .npz file for the Reso spectrometer
13
14    Args:
15        calibration_data_filepath (StrOrBytesPath): path of the file to load
16    """
17    try:
18        calibration = np.load(calibration_data_filepath)
19        wavelength_calibration = calibration["wavelength"]
20        lamp_spec = calibration["lamp_ref"]
21        lamp_measbyReso = calibration["lamp_measured"]
22        calibration_smoothed = np.abs(calibration["corr_factor_smoothed"])
23    except FileNotFoundError:
24        print(
25            "Error: calibration data for the UV spectrometer not found.\n"
26            "Please copy the folder Attoworld/src/attoworld/spectrum/calibration_data into your current working directory\n"
27            "or alternatively create a calibration file with relative path ./calibration_data/Reso_Spectrometer_CalibrationCorrection.npz\n"
28        )
29        raise FileNotFoundError("calibration data not found")
30    return wavelength_calibration, lamp_spec, lamp_measbyReso, calibration_smoothed

165def read_spectrometer_excel(filename):
166    """Reads xls file (passed as the string filename without extension!) produced by the UV-VIS spectrometer RESONANCE VS-7550-VUV.
167
168    This function currently only works in Linux (Ubuntu), since it uses the OS system command to brutally copy the content of the xls file into a txt;
169    such copying is necessary because the xls file is not readable by pandas (it's a non-standard xls, I couldn't find any other workaround).
170    For windows there is a similar command, please replace the line
171    os.system("cat " + filename + ".xls > " + filename + ".txt")
172
173    44 (skipped) rows at the beginning of the file (all headers)
174
175    Args:
176        filename: name of the file (without extension)
177
178    Returns:
179        numpy array with the spectral data (see original excel for more details)"""
180    if "." in filename and "./" not in filename:
181        raise ValueError(
182            "in function read_spectrometer_excel, filename must be passed without extension (filename should not contain a dot)"
183        )
184    os.system("cat " + filename + ".xls > " + filename + ".txt")
185    dataF = pandas.read_table(
186        filename + ".txt", sep="\t", keep_default_na=False, skiprows=44
187    )  # keep_default_na=False,
188    data = []
189    for row in dataF.values:
190        data.append([])
191        for x in row:
192            if x == "":
193                data[-1].append(float("nan"))
194            else:
195                data[-1].append(float(x))
196    return np.array(data)

 92def read_spectrum_ocean_optics(filename):
 93    """Reads a spectrum file from the Ocean Optics spectrometer and returns the wavelength and spectrum arrays.
 94    The file is expected to have one column of wavelength data and one column of spectrum data.
 95    wavelengths are in nm.
 96    14 lines of header are skipped."""
 97    data = pandas.read_table(filename, sep="\t", keep_default_na=True, skiprows=14)
 98    wvl = np.array(data.iloc[:, 0])
 99    spectrum = np.array(data.iloc[:, 1])
100    return wvl, spectrum

199def calibrate(
200    data,
201    column: int,
202    calibration_file_path="./calibration_data/Reso_Spectrometer_CalibrationCorrection.npz",
203    dark=None,
204    dark_c=None,
205    stitch: bool = False,
206    smooth_points: int = 10,
207    null_calibration: bool = False,
208    wavelength_calibration_intercept: float = 3.538,
209    wavelength_calibration_slope: float = 1.003427,
210):
211    """For the UV spectrometer. Calibrate the spectrum number 'column' in the the 'data' array. Intensity calibration factor is loaded in the first lines of UVSpectrumAnalysis.py.
212        Notice: each saved spectrum has 7 columns in the excel file; however, the argument column refers to the index of the spectrum, not the actual column in the file.
213
214    Args:
215        data: data [numpy] array (see read_spectrometer_excel)
216        column: index of the spectrum to be calibrated
217        calibration_file_path: path of the calibration file
218        dark: dark spectrum (optional)
219        dark_c: column of the dark spectrum (optional)
220        stitch: Boolean, if True, the spectra are stitched together (optional, default False)
221                if True, the spectrum will be calculated by stitching the 5 partial spectra in the previous columns.
222                the stitching only works if the spectrum number 'column' is a FULL-SPECTRUM in the xls file.
223                Stitching was implemented because the internal software of the spectrometer is not good enough at it (the saved full spectrum typically has discontinuities).
224        smooth_points: number of points for smoothing (optional, default 10). if 0, no smoothing is applied
225        null_calibration: Boolean, if True, the intensity calibration is not applied (optional, default False)
226        wavelength_calibration_intercept: intercept of the wavelength calibration (optional, default 3.538). If None no wavelength calibration is applied
227        wavelength_calibration_slope: slope of the wavelength calibration (optional, default 1.003427). If None no wavelength calibration is applied
228                Wavelength calibration is in the form λ_true = λ_measured * slope + intercept
229                Notice: wavelength calibration has an error of +/- 1  to 2 nm due to the stitching process.
230                    for a better calibration compute and apply the correction to each of the 5 partial spectra
231
232    Returns:
233        wavelength: wavelength array
234        spectrum: calibrated spectrum
235    """
236    wavelength_calibration, lamp_spec, lamp_measbyReso, calibration_smoothed = (
237        load_calibration_data(calibration_file_path)
238    )
239
240    if wavelength_calibration_slope is None or wavelength_calibration_intercept is None:
241        cal_slope = 1.0
242        cal_intercept = 0.0
243    else:
244        cal_slope = wavelength_calibration_slope
245        cal_intercept = wavelength_calibration_intercept
246
247    if stitch:
248        # check that the chosen spectrum (in column) is actually a full one
249        theoretically_complete_wvl = (
250            data[:, 7 * (column) + 1] * cal_slope + cal_intercept
251        )
252        theoretically_complete_wvl = theoretically_complete_wvl[
253            ~np.isnan(theoretically_complete_wvl)
254        ]
255        if theoretically_complete_wvl[-1] - theoretically_complete_wvl[0] < 900:
256            raise ValueError(
257                "in function calibrate(), the spectrum to be stitched is not a full spectrum, please check the column index"
258            )
259
260        extrema = []
261        spectra = []
262        for i in range(5):
263            wvl_single_sp = (
264                data[:, 7 * (column - 5 + i) + 1] * cal_slope + cal_intercept
265            )
266            single_sp = data[:, 7 * (column - 5 + i) + 2]
267            if dark is not None and dark_c is not None:
268                single_sp -= dark[:, 7 * (dark_c - 5 + i) + 2]
269            wvl_single_sp = wvl_single_sp[~np.isnan(wvl_single_sp)]
270            single_sp = single_sp[~np.isnan(single_sp)]
271            extrema.append([wvl_single_sp[0], wvl_single_sp[-1]])
272            spectra.append([wvl_single_sp, single_sp])
273
274        overlaps = []
275        centers_of_overlap = [extrema[0][0]]
276        for i in range(4):
277            overlaps.append(extrema[i][1] - extrema[i + 1][0])
278            centers_of_overlap.append(extrema[i][1] / 2 + extrema[i + 1][0] / 2)
279
280        w = np.min(overlaps)
281        centers_of_overlap.append(extrema[4][1] + w / 2)
282        centers_of_overlap[0] = centers_of_overlap[0] - w / 2
283
284        wavelength = data[:, 7 * column + 1] * cal_slope + cal_intercept
285        spectrum = wavelength * 0.0
286
287        for i in range(5):
288            center = centers_of_overlap[i] / 2 + centers_of_overlap[i + 1] / 2
289            FWHM = centers_of_overlap[i + 1] - centers_of_overlap[i]
290            spectrum += np.interp(
291                wavelength,
292                spectra[i][0],
293                spectra[i][1] * tukey_window(spectra[i][0], center, FWHM, w),
294            )
295            plt.plot(
296                spectra[i][0],
297                spectra[i][1] * tukey_window(spectra[i][0], center, FWHM, w),
298            )
299        calibr_interpolated = np.interp(
300            wavelength, wavelength_calibration, calibration_smoothed
301        )
302        if not null_calibration:
303            spectrum *= calibr_interpolated
304
305        plt.plot(
306            data[:, 7 * column + 1] * cal_slope + cal_intercept, data[:, 7 * column + 2]
307        )
308        plt.plot(wavelength, spectrum)
309
310        if w < 0:
311            raise ValueError(
312                "in function calibrate(), the spectra to be stitched do not overlap"
313            )
314
315    else:
316        wavelength = data[:, 7 * column + 1] * cal_slope + cal_intercept
317        if dark is None or dark_c is None:
318            calibr_interpolated = np.interp(
319                wavelength, wavelength_calibration, calibration_smoothed
320            )
321            if null_calibration:
322                calibr_interpolated = calibr_interpolated * 0.0 + 1
323            spectrum = (
324                data[:, 7 * column + 2] - np.min(data[:, 7 * column + 2])
325            ) * calibr_interpolated
326        else:
327            calibr_interpolated = np.interp(
328                wavelength, wavelength_calibration, calibration_smoothed
329            )
330            if null_calibration:
331                calibr_interpolated = calibr_interpolated * 0.0 + 1
332            spectrum = (
333                data[:, 7 * column + 2] - dark[:, 7 * dark_c + 2]
334            ) * calibr_interpolated
335
336            if len(data[:, 7 * column + 2][~np.isnan(data[:, 7 * column + 2])]) != len(
337                dark[:, 7 * dark_c + 2][~np.isnan(dark[:, 7 * dark_c + 2])]
338            ):
339                raise ValueError(
340                    "In function calibrate(), the dark spectrum is not the same length as the spectrum to be calibrated.\n"
341                    "Please check the column indices of the dark spectrum and the main spectrum"
342                )
343            if (
344                data[:, 7 * column + 1][~np.isnan(data[:, 7 * column + 1])][0]
345                != dark[:, 7 * dark_c + 2][~np.isnan(dark[:, 7 * dark_c + 1])][0]
346            ):
347                print(
348                    "WARNING: In function calibrate(), the dark spectrum does not have the same wavelength as the spectrum to be calibrated.\n"
349                    "Please check the column indices of the dark spectrum and the main spectrum"
350                )
351
352    if smooth_points is not None and smooth_points > 0:
353        spectrum = smooth(spectrum, smooth_points)
354    return wavelength, spectrum