Dear Maria, are you using a particular package for the deconvolution?

I suggest to check scipy deconvolution.

You should load your ascii file into a numpy array, for instance by using loadtxt. And define the shape of your filter.

Your question is too vague as for me to try to answer in a better fashion.

Sorry for being unspecific. Yes, I use loadtxt to make two arrays - for transmission, and for wavenumbers. And using pyplot I can create the plot, however it looks quite jagged and I read that Gaussian or Lorenzian deconvolution should be applied to make it more tidy. I assumed that I should use the transmission values as the x in a Gaussian distribution.

I suggest you to check the meaning of deconvolution.

The idea is that when you  measure, you convolute "the physics" with your instrument's resolution. So, something that is infinitely narrow in your spectra may appear as having a gaussian shape.

Therefore the first step you have to achieve is knowing how your instrument behaves. Take the spectra of something as narrow as you can get (a laser line, for instance) and use that measurement in your deconvolution function.

The help in scipy is here:

https://github.com/scipy/scipy/blob/v0.16.0rc1/scipy/signal/signaltools.py

(Newer commits than  the master branch have slightly better docstrings).

"""Deconvolves ``divisor`` out of ``signal``.

Returns the quotient and remainder such that

``signal = convolve(divisor, quotient) + remainder``

Parameters

----------

signal : array_like

Signal data, typically a recorded signal

divisor : array_like

Divisor data, typically an impulse response or filter that was

applied to the original signal

Returns

-------

quotient : ndarray

Quotient, typically the recovered original signal

remainder : ndarray

Remainder

Examples

--------

Deconvolve a signal that's been filtered:

>>> from scipy import signal

>>> original = [0, 1, 0, 0, 1, 1, 0, 0]

>>> impulse_response = [2, 1]

>>> recorded = signal.convolve(impulse_response, original)

>>> recorded

array([0, 2, 1, 0, 2, 3, 1, 0, 0])

>>> recovered, remainder = signal.deconvolve(recorded, impulse_response)

>>> recovered

array([ 0., 1., 0., 0., 1., 1., 0., 0.])

See also

--------

numpy.polydiv : performs polynomial division (same operation, but

also accepts poly1d objects)

"""

Perhaps Maria is only looking for a smooth interpolation for the data.

I'm assuming this because you've mentioned data comes out jagged using (matplotlib?) pyplot. Most plot libraries use a straight line between any two points, so if you lack data (step on your function may be too big) the plot will come out with that effect. By decreasing the step (which means increasing the number of samples, typically by interpolating between them) you can solve that issue, specially if you use methods adequate to do so. If this is the case I would advise Scipy "interp1d":

http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.interpolate.interp1d.html

http://docs.scipy.org/doc/scipy/reference/tutorial/interpolate.html

, perhaps using kind='cubic' interpolation. Please notice that interpolated data should only be assumed physically relevant if there's a reason for it. For plotting this is OK, for experimental procedures this is considered an increase in uncertainty (and rightfully so).

I'm not sure in what scientific area (and in what procedure) you're retrieving data from but it seems you're assuming the raw data already comes convoluted (which I find a bit unlikely). Convolution "per se" is a mathematical operation between two functions and you basically do it when you want to transform function 1 into function 3 by using function 2 (so f1*f2 = f3 ; actually I believe this is regardless of obtaining gaussian shape or not Aquiles, that will depend on the purpose).

I'm only saying this because by doing deconvolution you are hopping to recover one of the functions (let's say function 1 by using function 2 and 3). That will only solve the "jagged" behavior if that more "primal" signal is less "jagged" than the one you are dealing now. Again it seems to me unlikely in most uses I've seen of convolution (so convoluted signal is usually less jagged than original). I've just checked and the wikipedia article about convolution has some pretty nice animations explaining illustrating the effect of this operation:

https://en.wikipedia.org/wiki/Convolution

Thank you both, and apologies for the late reply. Indeed, it is interpolation that I needed but it seems that I cannot get it right within Python. In the interpolated graph I end up with values of -60 000 for intensity. But anyway, it was helpful to read your advice.

I was trying to make a better-looking plot for the theoretical vibration spectrum and according to Gabedit, I can pick Gaussian convolution or Lorentzian convolution to make it look like a proper spectrum instead of ugly straight lines between points.

I'm attaching an example. Hope it is clear now. 

You can also play around changing the value in the first linspace; if you go from 20 to, for instance 10, the effect is going to be interesting ;-)