#include <math.h>
#include <cxcore.h>
#include <cv.h>
#include <highgui.h>
#pragma comment(lib,"cv.lib")
#pragma comment(lib,"cxcore.lib")
#pragma comment(lib,"highgui.lib")

void edgeSmall(unsigned char *out, unsigned char *in, int inHeight, int inWidth, int inChannel, unsigned char *vec);
void sobelVec(unsigned char *out, unsigned char *vec, unsigned char *in, int inHeight, int inWidth, int inChannel);

int main(){
    IplImage *in, *out;
    unsigned char *vec;

    in = cvLoadImage("./image.bmp", 0);
    out = cvCreateImage(cvSize(in->width, in->height), in->depth, in->nChannels);
    vec = (unsigned char*)malloc(sizeof(unsigned char) * in->widthStep * in->height);
    
    sobelVec((unsigned char*)out->imageData, vec, (unsigned char*)in->imageData, in->height, in->width, in->nChannels);
    edgeSmall((unsigned char*)out->imageData, (unsigned char*)in->imageData, in->height, in->width, in->nChannels, vec);

    cvNamedWindow("in", 1);
    cvShowImage("in", in);
    cvNamedWindow("out", 1);
    cvShowImage("out", out);
    cvWaitKey(-1);
    
    cvDestroyWindow("in");
    cvDestroyWindow("out");

    cvReleaseImage( &in );
    cvReleaseImage( &out );
    
    free(vec);

    return 0;
}

void edgeSmall(unsigned char *out, unsigned char *in, int inHeight, int inWidth, int inChannel, unsigned char *vec){
    int i, j;
    int index;
    
    int widthStep;
    widthStep = inWidth * inChannel;
    if(widthStep % 4 != 0){
        widthStep = widthStep + 4 - widthStep % 4;
    }

    for(i = 1; i < (inHeight - 1); ++i){
        for(j = 1; j < (inWidth - 1); ++j){
            index = i * widthStep + j;
            
            switch(vec[index]){
            case 0:
                // 0度方向
                if(in[index] > in[index - 1] &&
                    in[index] > in[index + 1]){
                        out[index] = in[index];
                } else {
                    out[index] = 0;
                }
                break;
            case 1:
                // 45度方向
                if(in[index] > in[index - 1 + widthStep] &&
                    in[index] > in[index + 1 - widthStep]){
                        // 45度方向
                        out[index] = in[index];
                } else {
                    out[index] = 0;
                }
                break;
            case 2:
                // 90度方向
                if(in[index] > in[index + widthStep] &&
                    in[index] > in[index - widthStep]){
                        out[index] = in[index];
                } else {
                    out[index] = 0;
                }
                break;
            case 3:
                // 135度方向
                if(in[index] > in[index + 1 + widthStep] &&
                    in[index] > in[index - 1 - widthStep]){
                        out[index] = in[index];
                } else {
                    out[index] = 0;
                }
                break;
            }
        }
    }
}

void sobelVec(unsigned char *out, unsigned char *vec, unsigned char *in, int inHeight, int inWidth, int inChannel){
    int i, j, k, l;
    int weightH[9] = {-1, 0, 1,
                     -2, 0, 2,
                     -1, 0, 1};
    int weightV[9] = {-1, -2, -1,
                      0,  0,  0,
                      1,  2,  1};
    int *imgTmp;
    int weightSize = 3;

    int start = -1;
    int end = start + weightSize;
    double dataTmpH;
    double dataTmpV;
    double dataTmpSum;
    double vecTan;
    int widthStep;

    widthStep = inWidth * inChannel;
    if(widthStep % 4 != 0){
        widthStep = widthStep + 4- widthStep%4;
    }

    for(i = 1; i < (inHeight - 1); ++i){
        for(j = 1; j < (inWidth - 1); ++j){
            dataTmpH = 0;
            dataTmpV = 0;
            for(k = start; k < end; ++k){
                for(l = start; l < end; ++l){
                     //水平方向
                    dataTmpH += weightH[(k - start) * weightSize + (l - start)] * (unsigned char)in[(i + k) * widthStep + (j + l)];

                     //垂直方向
                    dataTmpV += weightV[(k - start) * weightSize + (l - start)] * (unsigned char)in[(i + k) * widthStep + (j + l)];

                }
            }
            dataTmpSum = sqrt(dataTmpH * dataTmpH + dataTmpV * dataTmpV);
            vecTan = dataTmpV / dataTmpH;

            if(-0.4142 < vecTan && vecTan <= 0.4142){
                // 0度方向
                vec[i * widthStep + j] = 0;
            } else if(0.4142 < vecTan && vecTan < 2.4142){
                // 45度方向
                vec[i * widthStep + j] = 1;
            } else if(abs(vecTan) >= 2.4142){
                // 90度方向
                vec[i * widthStep + j] = 2;
            } else{
                // 135度方向
                vec[i * widthStep + j] = 3;
            }
            
            //255を超えた値は255、0未満の値は絶対値
            if(dataTmpSum > 255){
                dataTmpSum = 255;
            } else if(dataTmpSum < 0){
                dataTmpSum = -dataTmpSum;
            }
            out[i * widthStep + j] = dataTmpSum;
        }
    }

}