Sinc-like feature removal with an FFT

Outline

I have some Raman spectroscopy data that contains features that I want to remove. I'll point them out in the figure below:

Here you can see two sets of peaks that I am interested in, with the right set annotated. The main "vibrational" peaks (blue) I want to preserve, and the accompanying "rotational" peaks (red) that I want to remove with some kind of filter. There is also a potential for other vibrational peaks to appear (green) among the rotational peaks that I also want to preserve.

What I've tried

My idea is that, because the rotational peaks appear to have a regular frequency akin to a sinc function, I can remove them by taking an FFT, removing the associated frequencies, and retrieving the filtered spectrum with an IFFT.

Here is the result of my attempt:

This figure shows a crop of the region annotated in the previous figure. Hopefully the figure legends are clear, but I will point out some key things:

• I first subtracted the minimum value for the cropped region, to make the filtering process simpler
• I have used a low-pass filter (LPF) to remove the Gaussian-like shape from the FFT (bottom; red) of the original spectrum (top; black) to produce the filtered FFT (bottom; green)
• I perform an IFFT on the filtered FFT to produce a spectrum containing my rotational peaks (top; green), and then subtract it from my original spectrum to obtain my filtered result (top; magenta).

You can see an overlaid comparison of the original and filtered spectra below:

Thoughts/Observations

As you can see my filtering attempt is not perfect, however the result of the LPF and IFFT does seem to do a reasonable job at recreating the rotational peaks and not the vibrational peaks, and the peak of interest (encircled green in the first figure) is made more prominent.

The main issues from what I can see are:

• The LPF has not perfectly fit to the Gaussian-like shape (although I am hesitant to fine-tune this value as then the process will lose some generalisability), this has caused the rotational peaks to not be fully removed.
• The filtered spectrum seems to have added a slight sine wave across the filtered regions, which I've roughly annotated
• The left filtered region is overall lower in intensity than before, which could be easily fixed with an additional step afterwards, but I think it could prove unnecessary with some improvements to the filtering process

Code & Data

I'll post a code snippet for how I have applied my filter:

import numpy as np
from scipy.fft import fft, fftshift, ifft, ifftshift

# Make a copy of the spectrum
spectrum_filtered = np.copy(spectrum_avg)
# Manually define some regions to be cropped and filtered
crop_regions = [[379, 479], [620, 791]]
for crop_region in crop_regions:
    # Crop the spectrum
    crop_left, crop_right = crop_region
    crop = np.copy(spectrum_avg[crop_left:crop_right])
    
    # Subtract the minimum value
    min_val = np.min(crop)
    crop -= min_val
    
    # FFT
    crop_fft = np.abs(fftshift(fft(crop)))
    
    # Fit the baseline using a low-pass filter
    fit = butter_filter(crop_fft, cutoff=3, fs=crop_fft.size, order=5, btype='lowpass')
    
    # Remove the baseline, then perform an IFFT (this holds the rotational peaks)
    crop_fft_fit = crop_fft - fit
    crop_ifft = np.abs(ifftshift(ifft(crop_fft_fit)))
    
    # Subtract the rotational peaks
    crop_filtered = crop - crop_ifft

    # Replace the cropped region with the filtered spectrum (re-adding the minimum value)
    spectrum_filtered[crop_left:crop_right] = crop_filtered + min_val

Here is a pastebin link to the example spectrum: https://pastebin.com/HHuJJXxL

I'm happy to clarify anything should it prove necessary. I would greatly appreciate any feedback, criticism, or suggestions on how to improve this method!

It caught my interest as a model fit problem so I processed the spectrum around the 600 to 900 region. Whatever the shape is it cannot be reasonably be described as a sinc function which is double width at the origin. The rotational modes are zero at the origin.

I think the correct empirical functional form to model this unwanted spectrum of rotational peaks around centre frequency channel X0 is roughly as follows (my quantum mechanics is quite rusty).

 f(x) = A+B*abs(x-X0)*exp(-(x-X0)^2/C)*sin((x-X0)*D)

For the main peak and only to the nearest channel I get maximum as2.018705692941782763e+04 in channel x0 = 707. The model parameters for the fit are approximate but are precise enough to show the asymmetry of the peak profile either side of the main peak.

BTW No attempt made at baseline removal.

The values are just a starting guide you would have to fit them properly with least squares to the local spectrum along with baseline removal. I have fitted to the wings away from the Raman zone. Not entirely happy with the high side fit but it is probably not far off.

The high side peak fit needs additional work. Is there some interference from Raman scattering near the main peak?