From f6f5e938c3d6949ed435512ae936f98b827604a4 Mon Sep 17 00:00:00 2001 From: marcin-witkowski <59771765+marcin-witkowski@users.noreply.github.com> Date: Thu, 9 Feb 2023 21:22:49 +0100 Subject: [PATCH] Component selection in Raman and spectrum saving added Added the ability to select the component of the polarizability tensor for calculations of the Raman cross section. Added the ability to save the generated spectrum in xy format. Fixed some typos and errors in the process. --- plot_RamanIR.py | 48 ++++++++++++++++++++++++++++++++++++------------ 1 file changed, 36 insertions(+), 12 deletions(-) diff --git a/plot_RamanIR.py b/plot_RamanIR.py index df2ee04..9663203 100755 --- a/plot_RamanIR.py +++ b/plot_RamanIR.py @@ -1,6 +1,6 @@ #! /usr/bin/env python from chemPackage import collect, dressedT -from numpy import linspace +from numpy import linspace, savetxt, array from chemPackage.constants import PI import argparse import matplotlib.pyplot as plt @@ -19,7 +19,7 @@ def main(): "molecule Raman scattering or DIM/QM SERS.") parser.add_argument("--dressed", "-d", nargs=2, metavar=("dim.out", "freq.out"), - type=str, help="Plot the raman spectrum using dressed tensors formulism. " + type=str, help="Plot the raman spectrum using dressed tensors formalism. " "Requires the output from a DIM calculation for the field, and the output " "of a frequency calculation. ", dest="dressed", default=None) @@ -37,6 +37,9 @@ def main(): parser.add_argument("--fwhm", metavar="20", default=20, type=float, help="Full width at half max for the peaks in the spectrum") + + parser.add_argument("--component", "-c", metavar="xx", default="all", type=str, + help="Supply the component of the tensor's derivative to calculate Raman intensity") parser.add_argument("--scalefactor", metavar="30", default=30, type=float, help="This is a scale factor for the y axis.") @@ -52,11 +55,12 @@ def main(): parser.add_argument("--figname", default=None, type=str, metavar="fig.png", help="Supply a name for the figure if it is to be saved") + parser.add_argument("--specname", default=None, type=str, metavar="spectrum.txt", + help="Supply a name for the spectrum file if it is to be saved in xy format") + parser.add_argument("--no-sticks", default=True, action='store_false',dest='lsticks', help="If included sticks will not be plotted") - - args = parser.parse_args() # sanity check: did they choose raman or dressed tensors. If not, give up on them @@ -92,7 +96,7 @@ def main(): # plot that spectrum dressedTensors(dimout, freqout, args.translation, args.fwhm, args.scalefactor, - args.xaxis, args.yaxis,args.lsticks, args.figname) + args.xaxis, args.yaxis,args.lsticks, args.figname, args.specname) @@ -104,7 +108,8 @@ def main(): exit("Please use the --help option for more information") # plot the spectrum - plotRaman(freqout, args.fwhm, args.scalefactor, args.xaxis, args.yaxis, args.lsticks, args.figname) + plotRaman(freqout, args.fwhm, args.scalefactor, args.xaxis, args.yaxis, args.component, + args.lsticks, args.figname, args.specname) # The second real condition: Raman by numerical differentiation. elif args.IR != None and args.dressed == None: @@ -114,7 +119,7 @@ def main(): exit("Please use the --help option for more information") # plot the spectrum - plotIR(freqout, args.fwhm, args.xaxis, args.yaxis, args.lsticks, args.figname) + plotIR(freqout, args.fwhm, args.xaxis, args.yaxis, args.lsticks, args.figname, args.specname) else: exit("This shouldn't even be possible, I guess I'm not as intellegent as I thought I was." @@ -124,7 +129,7 @@ def main(): -def dressedTensors(dimout, freqout, tr, fwhm, scaleexp, xaxis, yaxis, lsticks, figname): +def dressedTensors(dimout, freqout, tr, fwhm, scaleexp, xaxis, yaxis, lsticks, figname, specname): """ This function plots the Raman spectrum of a system using the Dressed Tensors formulism. @@ -143,6 +148,9 @@ def dressedTensors(dimout, freqout, tr, fwhm, scaleexp, xaxis, yaxis, lsticks, f the form of a list. ie. [lower, upper] figname ==> If it isn't None then this function will save the figure to filename supplied. + specname ==> If it isn't None then this function will save the spectrum in the + xy format to the filename supplied. + """ # collect the dipole-dipole polarizabilities derivatives and @@ -181,15 +189,18 @@ def dressedTensors(dimout, freqout, tr, fwhm, scaleexp, xaxis, yaxis, lsticks, f # check if ranges for the specturm were given. If given apply them sub.set_xlim(xaxis[0], xaxis[1]) if yaxis != None: - sub.set_ylim(yaxix[0], yaxis[1]) + sub.set_ylim(yaxis[0], yaxis[1]) if figname != None: plt.savefig(figname, dpi=300) else: plt.show() + + if specname != None: + savetxt(specname, array([domain, y2*scale]).T, fmt=['%.2f', '%4.8g']) -def plotRaman(freqout, fwhm, scaleexp, xaxis, yaxis, lsticks, figname): +def plotRaman(freqout, fwhm, scaleexp, xaxis, yaxis, component, lsticks, figname, specname): """ This function plots the Raman spectrum of a system using numerical differentiation around the normal modes. @@ -202,8 +213,12 @@ def plotRaman(freqout, fwhm, scaleexp, xaxis, yaxis, lsticks, figname): the form of a list. ie. [lower, upper] yaxis ==> gives the lower and upper bounds for the y-axis. It is supplied in the form of a list. ie. [lower, upper] + component ==> gives the component of the polarizability tensor for which the + Raman spectrum is to be calculated figname ==> If it isn't None then this function will save the figure to filename supplied. + specname ==> If it isn't None then this function will save the spectrum in the + xy format to the filename supplied. """ # good news, because the chem package is so great this works for both # the Raman of the free molecule and SERS calculated with DIM/QM @@ -211,7 +226,7 @@ def plotRaman(freqout, fwhm, scaleexp, xaxis, yaxis, lsticks, figname): # collect the dipole-dipole polarizabilities derivatives freqout.collect_raman_derivatives() # calculate the raman cross section for the system. - raman_intensity = freqout.cross_section() + raman_intensity = freqout.cross_section(component=component) scale = 10**scaleexp @@ -241,8 +256,11 @@ def plotRaman(freqout, fwhm, scaleexp, xaxis, yaxis, lsticks, figname): plt.savefig(figname, dpi=300) else: plt.show() + + if specname != None: + savetxt(specname, array([domain, y2*scale]).T, fmt=['%.2f', '%4.8g']) -def plotIR(freqout, fwhm, xaxis, yaxis, lsticks, figname): +def plotIR(freqout, fwhm, xaxis, yaxis, lsticks, figname, specname): """ This function plots the Raman spectrum of a system using numerical differentiation around the normal modes. @@ -253,6 +271,8 @@ def plotIR(freqout, fwhm, xaxis, yaxis, lsticks, figname): the form of a list. ie. [lower, upper] figname ==> If it isn't None then this function will save the figure to filename supplied. + specname ==> If it isn't None then this function will save the spectrum in the + xy format to the filename supplied. """ fig = plt.figure() @@ -278,6 +298,9 @@ def plotIR(freqout, fwhm, xaxis, yaxis, lsticks, figname): else: plt.show() + if specname != None: + savetxt(specname, array([domain, y]).T, fmt=['%.2f', '%4.8g']) + def lorentzian(x, peak=0, height=1.0, fwhm=None, hwhm=None): '''Calculates a three-parameter lorentzian for a given domain.''' @@ -312,6 +335,7 @@ def sum_lorentzian(x, peak=None, height=None, fwhm=None, hwhm=None): try: main() except KeyboardInterrupt: + import sys sys.exit(1)