diff --git a/mcutils/mcutils.py b/mcutils/mcutils.py
index 58db3fc..0c10db7 100644
--- a/mcutils/mcutils.py
+++ b/mcutils/mcutils.py
@@ -288,7 +288,7 @@ def convolveQuad(x,y,xres,yres,useAutoCrop=True,
     x_integral = xresInv+x[i]
     ytemp = interpolate(x,y,x_integral,fill_value=fill_value,axis=axis)/area
     # do integration
-    res = integrate.simps(ytemp*yresInv,x=xresInv+x[i],axis=axis)
+    res = integrate.simps(ytemp[0]*yresInv,x=xresInv+x[i],axis=axis)
     colon = (slice(None),)
     setElement(out,i,res,axis=axis)
   return out
@@ -326,6 +326,14 @@ def convolve(x,y,xres,yres,useAutoCrop=True,approximantOrder=None,fill_value=0.,
   
   ## prepare output
   out = np.empty_like(y)
+
+  ## fill up NaN
+  for i in range(y.shape[0]):
+    isok = np.isfinite(y[i])
+    xp   = x[isok]
+    fp   = y[i,isok]
+    y[i] = np.interp(x, xp, fp)
+
   nP = len(x)
   f = interpolator(x,y,fill_value=fill_value,axis=axis)
   for i in range(nP):
@@ -335,10 +343,10 @@ def convolve(x,y,xres,yres,useAutoCrop=True,approximantOrder=None,fill_value=0.,
     #ytemp  = interpolate(x,y,x_integral,fill_value=fill_value,axis=axis)/area
     #ytemp = f(x_integral)/area
     # do integration
-    res = integrate.simps( np.transpose(ytemp.T*yresInv),x=xresInv+x[i],axis=axis)
+    res = integrate.simps( ytemp*yresInv,x=xresInv+x[i],axis=axis)
     colon = (slice(None),)
     setElement(out,i,res,axis=axis)
-  return out
+  return out.T
 
 
 def fftconvolve(x,y,yres,xres=None,normalize=False):