How to remove the boundary effects arising due to zero padding in scipy/numpy fft? I have made a python code to smoothen a given signal using the Weierstrass transform, which is basically the convolution of a normalised gaussian with a signal.

The code is as follows:  



#Importing relevant libraries  
from __future__ import division  
from scipy.signal import fftconvolve   
import numpy as np 

def smooth_func(sig, x, t= 0.002):   
    '''
    x is an 1-D array, 
    sig is the input signal and a function of x.
    '''
    N = len(x)   
    x1 = x[-1]   
    x0 = x[0]    


# defining a new array y which is symmetric around zero, to make the gaussian symmetric.  
    y = np.linspace(-(x1-x0)/2, (x1-x0)/2, N)   
    #gaussian centered around zero.  
    gaus = np.exp(-y**(2)/t)     

#using fftconvolve to speed up the convolution; gaus.sum() is the normalization constant.  
    return fftconvolve(sig, gaus/gaus.sum(), mode='same')




If I run this code for say a step function, it smoothens the corner, but at the boundary it interprets another corner and smoothens that too, as a result giving unnecessary behaviour at the boundary. I explain this with a figure shown in the link below.
Boundary effects

This problem does not arise if we directly integrate to find convolution. Hence the problem is not in Weierstrass transform, and hence the problem is in the fftconvolve function of scipy. 

To understand why this problem arises we first need to understand the working of fftconvolve in scipy.
The fftconvolve function basically uses the convolution theorem to speed up the computation.
In short it says:
convolution(int1,int2)=ifft(fft(int1)*fft(int2))
If we directly apply this theorem we dont get the desired result. To get the desired result we need to take the fft on a array double the size of max(int1,int2). But this leads to the undesired boundary effects. This is because in the fft code, if size(int) is greater than the size(over which to take fft) it zero pads the input and then takes the fft. This zero padding is exactly what is responsible for the undesired boundary effects.  

Can you suggest a way to remove this boundary effects?

I have tried to remove it by a simple trick. After smoothening the function I am compairing the value of the smoothened signal with the original signal near the boundaries and if they dont match I replace the value of the smoothened func with the input signal at that point.
It is as follows:  



i = 0 
eps=1e-3
while abs(smooth[i]-sig[i])&gt; eps: #compairing the signals on the left boundary
    smooth[i] = sig[i]
    i = i + 1
j = -1

while abs(smooth[j]-sig[j])&gt; eps: # compairing on the right boundary.
    smooth[j] = sig[j]
    j = j - 1




There is a problem with this method, because of using an epsilon there are small jumps in the smoothened function, as shown below:
jumps in the smooth func 

Can there be any changes made in the above method to solve this boundary problem?

Also I tried removing the zero padding in the fft source code and replaced it with a constant value, but it gave undesired results.
Can you suggest a way of removing this zero padding in the scipy fft source code?
Can there be any changes made in the above method to solve this
  boundary problem?


Yes there can be, but it all depends what you want to occur at the boundary. If you don't want it to be smoothed then simply replace the edge with the original signal for a sensible number of samples, depending on the width of your gaussian. You gaussian might be > 0 for say 10 samples, so replace the 5 samples at the start and end of the signal with the original samples.


  Can you suggest a way of removing this zero padding in the scipy fft
  source code?


Yes. Pad it yourself. If you want to consider the wrap around, as if your signal is truly periodic, then fabricate the wraparound effect for a power of 2 in size by carefully padding your own signal. Just as a simple idea, multiply your window by 2 and find the nearest power of 2, then fill the middle of the copy with zeros to pad it out

[1234567890] [12345] - zero pad here - [67890]


Now when you apply the convolution, you should get the same effect for the first window's worth of the result as if the convolution kernel had wrapped around the original signal. To illustrate, take your gaussian and zero pad it at the end, and slide through the convolution to see for yourself that it yields the same thing

 original   left    padding    right
[1234567890 12345 000000000000 67890]
[xgaussianx 00000 000000000000 00000]


so the value for the right most sample comes from this:

 original   left    padding    right
[1234567890 12345 000000000000 67890]
[sianx00000 00000 000000000000 xgaus]


Which is equivalent to what we actually want

[1234567890]
[xsianxgaus] (gaussian shifted by 5 samples and wrapped around)


and the value for the left most sample comes from

 original   left    padding    right
[1234567890 12345 000000000000 67890]
[00000xgaus sianx 000000000000 00000]


Which is also equivalent to what we actually want

[1234567890]
[sianxxgaus] (gaussian shifted by -5 samples and wrapped around)


All this relies on deciding what it is you actually want to happen. If you understand the mathematical process of convolution, and you understand exactly what it is you want to calculate, you can fabricate the conditions under power of 2 scenarios to ensure you get fast computation. Generally speaking most people just opt for the replace because their signal windows are not periodic and have sharp discontinuities at the edges.