Deriche edge detector

edit

I have read @Miki answer in this post and disable CV_OCL_RUN and CV_IPP_RUN in my own copy of opencv canny. results is here :

#include<opencv2/opencv.hpp>
#include "opencv2/core/ocl.hpp"
#include<iostream>
using namespace cv;


using namespace cv;

/*
Using Canny's Criteria to Derive a Recursively Implemented Optimal Edge Detector International Journal of Computer Vision,167-187 (1987)
http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.476.5736&rep=rep1&type=pdf
http://www.esiee.fr/~coupriem/Algo/algoima.html
*/


class ParallelGradientDericheYCols: public ParallelLoopBody
{
private:
    Mat &img;
    Mat &im1;
    double alphaMoyenne;
    double alphaDerive;
    bool verbose;

public:
    ParallelGradientDericheYCols(Mat& imgSrc, Mat &d,double alm,double ald):
        img(imgSrc),
        im1(d),
        alphaMoyenne(alm),
        alphaDerive(ald),
        verbose(false)
    {}
    void Verbose(bool b){verbose=b;}
    virtual void operator()(const Range& range) const
    {

        if (verbose)
            std::cout << getThreadNum()<<"# :Start from row " << range.start << " to "  << range.end-1<<" ("<<range.end-range.start<<" loops)" << std::endl;

        float               *f1,*f2;

        long                i,j,nb;
        long                tailleSequence=(img.rows>img.cols)?img.rows:img.cols;
        double              *g1=new double[tailleSequence],*g2=new double[tailleSequence];
        double  p=alphaDerive;
        double  k=pow(1-exp(-p),2.0)/(1+2*p*exp(-p)-exp(-2*p));
        double  kp=pow(1-exp(-p),2.0)/exp(-p);
        double a1=1,a2=0;
        double a3=0,a4=0;
        double b1=0,b3=2*exp(-p);
        double b2=0,b4=-2*exp(-p);
        int rows=img.rows,cols=img.cols;


        kp=pow(1-exp(-p),2.0)/exp(-p);
        a1=0,a2=kp*exp(-p),a3=-kp*exp(-p),a4=0;
        b1=2*exp(-p);
        b2=-exp(-2*p);


        switch(img.depth()){
        case CV_8U :
            {
            unsigned char       *c1;
            for (nb=0;nb<img.channels();nb++)
                {
                for (j=range.start;j<range.end;j++)
                    {
                    // Filtre causal vertical
                    c1 = (unsigned char*)img.ptr(0)+nb;
                    f2 = (float*)im1.ptr(0)+nb;
                    f2 += j;
                    c1+=j;
                    i=0;
                    g1[i] = (a1 + a2 +b1+b2)* *c1  ;
                    g1[i] = (a1 + a2 )* *c1  ;
                    i++;
                    c1+=cols;
                    g1[i] = a1 * *c1 + a2 * c1[-cols] + (b1+b2) * g1[i-1];
                    g1[i] = a1 * *c1 + a2 * c1[-cols] + (b1) * g1[i-1];
                    i++;
                    c1+=cols;
                    for (i=2;i<rows;i++,c1+=cols)
                        g1[i] = a1 * *c1 + a2 * c1[-cols] +b1*g1[i-1]+
                                b2 *g1[i-2];

                    // Filtre anticausal vertical
                    c1 = (unsigned char*)img.ptr(0)+nb;
                    c1 += (rows-1)*cols+j;
                    i = rows-1;
                    g2[i] =(a3+a4+b1+b2)* *c1; 
                    g2[i] =(a3+a4)* *c1; 
                    i--;
                    c1-=cols;
                    g2[i] = a3* c1[cols] + a4 * c1[cols]+(b1+b2)*g2[i+1];
                    g2[i] = a3* c1[cols] + a4 * c1[cols]+(b1)*g2[i+1];
                    i--;
                    c1-=cols;
                    for (i=rows-3;i>=0;i--,c1-=cols)
                        g2[i] = a3*c1[cols] +a4* c1[2*cols]+
                                b1*g2[i+1]+b2*g2[i+2];
                    for (i=0;i<rows;i++,f2+=cols)
                        *f2 = (float)(g1[i]+g2[i]);
                    }
                }
            }
            break;
        case CV_16S :
        case CV_16U :
            {
            unsigned short      *c1;
            for (nb=0;nb<img.channels();nb++)
                {
                for (j=range.start;j<range.end;j++)
                    {
                    c1 = ((unsigned short*)img.ptr(0))+nb;
                    f2 = ((float*)im1.ptr ...