How to export a calculated Raman spectrum in Q-Chem, as a Gaussian profile?

Question

Does anyone know how to export a calculated Raman spectrum (my calculation was performed in Q-Chem) as a data-file, which can be used to plot the profile of a Gaussian-shape spectrum (in origin, excel, etc.), rathe than centroid data points?

Answer 1

I wrote the following simple Python function that allows one to use any distribution from scipy.stats for peak broadening (the scale parameter determines how broad peaks will be):

import numpy as np
from scipy.stats import cauchy, norm


def broaden_spectrum(x, x0, y0, distribution="gaussian", scale=1.0, 
                     fit_points=True, *args, **kwargs):
    if distribution in {"gaussian", "norm"}:
        distribution = norm
    elif distribution in {"lorentzian", "cauchy"}:
        distribution = cauchy

    s = np.sum([yp * distribution.pdf(x, xp, scale=scale, *args, **kwargs)
                for xp, yp in zip(x0, y0)], axis=0)

    if fit_points:
        s_max = np.max(s)
        if s_max == 0.0:
            s_max = 1.0
        return s * np.max(y0) / s_max
    return s

As Geoff Hutchison mentioned, you can use cclib to read your logfile:

import cclib


data = cclib.ccopen("data/2-propanol.out").parse()
x = np.linspace(data.vibfreqs.min() - 100.,
                data.vibfreqs.max() + 100., num=1000)
y = broaden_spectrum(x, data.vibfreqs, data.vibirs, scale=40.)

Plotting is easy:

import matplotlib.pyplot as plt


plt.plot(x, y)
plt.gca().invert_xaxis()
plt.gca().invert_yaxis()
plt.ylabel("Absorbance (arb. unit.)")
plt.xlabel("Wavenumber (cm$^{-1}$)")
plt.title("Predicted IR spectrum of 2-propanol at revPBE-D4-gCP/def2-SVP")
plt.savefig("ir-spectrum.png")

This is easily customizable and the broaden_spectrum function works for all kinds of spectra.

Answer 2

Many thanks to Felipe, the code is very useful! I know this is just a simple editing, but this is how I adjusted part of the code for Raman spectra:

import cclib


data = cclib.ccopen("molecule.out").parse()
x = np.linspace(data.vibfreqs.min() - 100.,
                data.vibfreqs.max() + 100., num=1000)
y = broaden_spectrum(x, data.vibfreqs, data.vibramans, scale=10.)

And also this will save a .csv which can be used for further data processing:

import matplotlib.pyplot as plt
import pandas as pd


plt.plot(x, y)
plt.ylabel("Intensity (a.u.)")
plt.xlabel("Raman Shift (cm$^{-1}$)")
plt.title("Predicted Raman spectrum of molecule at B3LYP/6-311++G**")
plt.savefig("Raman_gauss.pdf", dpi=300, bbox_inches='tight')


dataframe = pd.DataFrame({'x': x,
                          'y':y})
dataframe.to_csv("Raman_spectrum.csv", index=False)

I'm sure this has multiple options for customization, but I though it might be useful.