Covolution function fitting problem

#pragma rtGlobals=3     // Use modern global access method and strict wave access.
Function Convolres(pw, yw, xw) : FitFunc    // linear DOS * Fermi, convolved with Gaussian energy resolution
    Wave pw, yw, xw
 
    // pw[0] = offset
    // pw[1] = energy resolution FWHM
    // pw[2] = square function initial position 
    // pw[3] = square function finish position 
    // pw[4] = square function amplitude
    
 
    
    // Make the resolution function wave W_res.
    Variable x0 = xw[0]
    Variable dx = (xw[inf]-xw[0])/(numpnts(xw)-1)                           // assumes even data spacing, which is necessary for the convolution anyway
    Make/FREE/D/N=(min(max(abs(3*pw[1]/dx), 5), numpnts(yw))) W_res // make the Gaussian resolution wave
    Redimension/N=(numpnts(W_res)+!mod(numpnts(W_res), 2)) W_res    // force W_res to have odd length
    SetScale/P x, -dx*(numpnts(W_res)-1)/2, dx, W_res
    W_res = gauss( x, 0, pw[1]/(2*sqrt(2*ln(2))) ) +pw[0]
    Variable a = sum(W_res)
    W_res /= a
 
    // To eliminate edge effects due to the convolution, add points to yw, make the spectrum,
    // do the convolution, and then delete the extra points.
    Redimension/N=(numpnts(yw)+2*numpnts(W_res)) xw, yw
    xw = (x0-numpnts(W_res)*dx) + p*dx
 
    // square function
    variable i
    yw=0+pw[0]
    for(i=pw[2];i<pw[3];i+=1)
    yw[i]=pw[4]
    endfor
    
    Convolve/A W_res, yw
    DeletePoints 0, numpnts(W_res), xw, yw
    DeletePoints numpnts(yw)-numpnts(W_res), numpnts(W_res), xw, yw
End