#include "stdafx.h"
#include "cv.h"
#include "cxcore.h"
#include "highgui.h"
#include "fft.h"
#define WINDOWNAME "Seam Carving"
#define GAUSSIAN_MASK_SIZE 7
#define SALIENCY_SIZE_W (img->width / 4)
#define SALIENCY_SIZE_H (img->height / 4)
#define BIAS_WEIGHT 10
#define threshold 0.01
#define iceil(expr) (int)(ceil((double)expr))
IplImage *f ,
*inv ,
*re ,
*im ,
*re2 ,
*im2 ,
*modular ,
*sine ,
*cosine ,
*logamp ,
*resi,
*tmp, *tmp2;
CvMat *sce, *sch16, *sch8, *sch4, *sch2;
void fft2(IplImage *src, IplImage *dst)
{ //实部、虚部
IplImage *image_Re = 0, *image_Im = 0, *Fourier = 0;
// int i, j;
image_Re = cvCreateImage(cvGetSize(src), src->depth, 1); //实部
//Imaginary part
image_Im = cvCreateImage(cvGetSize(src), src->depth, 1); //虚部
//2 channels (image_Re, image_Im)
Fourier = cvCreateImage(cvGetSize(src), src->depth, 2);
//实部的值初始设为源图像，虚部的值初始设为0
// Real part conversion from u8 to 64f (double)
cvConvertScale(src, image_Re, 1, 0);
// Imaginary part (zeros)
cvZero(image_Im);
// Join real and imaginary parts and stock them in Fourier image
cvMerge(image_Re, image_Im, 0, 0, Fourier);
// Application of the forward Fourier transform
cvDFT(Fourier, dst, CV_DXT_FORWARD);
cvReleaseImage(ℑ_Re);
cvReleaseImage(ℑ_Im);
cvReleaseImage(&Fourier);
}
void init(IplImage *src, IplImage *src0) {
tmp = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
f = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 2),
inv = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 2),
re = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
im = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
re2 = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
im2 = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
modular = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
sine = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
cosine = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
logamp = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1),
resi = cvCreateImage(cvGetSize(src), IPL_DEPTH_32F, 1);
sce = cvCreateMat(iceil(src0->height / 16) * 16, iceil(src0->width / 16) * 16, CV_32F);
sch16 = cvCreateMat(iceil(src0->height / 16), iceil(src0->width / 16), CV_32F);
sch8 = cvCreateMat(iceil(src0->height / 8 ), iceil(src0->width / 8 ), CV_32F);
sch4 = cvCreateMat(iceil(src0->height / 4 ), iceil(src0->width / 4 ), CV_32F);
sch2 = cvCreateMat(iceil(src0->height / 2 ), iceil(src0->width / 2 ), CV_32F);
};
void spect_ampl(IplImage *src, IplImage* &salmap, double maxScale = 1.0)
{
cvConvertScale(src, tmp, 1.0/255.0, 0);
fft2(tmp, f);
cvSplit(f, re, im, NULL, NULL);
cvPow(re, re2, 2.0); cvPow(im, im2, 2.0); cvAdd(re2, im2, modular); cvPow(modular, modular, 0.5);
cvDiv(im, modular, sine); cvDiv(re, modular, cosine);
cvLog(modular, logamp);
cvSmooth(logamp, resi, CV_BLUR, 3, 3);
cvSub(logamp, resi, resi);
cvExp(resi, resi);
cvMul(resi, cosine, re2); cvMul(resi, sine, im2); cvMerge(re2, im2, NULL, NULL, f);
cvDFT(f, inv, CV_DXT_INV_SCALE);
cvSplit(inv, re, im, NULL, NULL);
cvPow(re, re2, 2.0); cvPow(im, im2, 2.0); cvAdd(re2, im2, modular);
cvSmooth(modular, modular, CV_GAUSSIAN, GAUSSIAN_MASK_SIZE, GAUSSIAN_MASK_SIZE);
// normalization
double min = 0, max = 0;
cvMinMaxLoc(modular, &min, &max, NULL, NULL);
max-=min;
cvConvertScale(modular, salmap, maxScale / max, -min / max * maxScale);
cvThreshold(salmap, salmap, threshold, threshold, CV_THRESH_TRUNC);
//cvConvertScale(salmap, salmap, 1.0 / threshold);
}
inline double sgetReal(CvMat *mat, int y, int x, double defval = 0) {
if (x < 0 || y < 0 || x >= mat->width || y >= mat->height) return defval;
return ((float*)(mat->data.ptr + (size_t)mat->step*y))[x];
}
inline void ssetReal(CvMat *mat, int y, int x, double val) {
if (x < 0 || y < 0 || x >= mat->width || y >= mat->height) return;
((float*)(mat->data.ptr + (size_t)mat->step*y))[x] = val;
}
inline double min3(double a, double b, double c) {
return (a < b) ? ((a < c) ? a : c) : ((b < c) ? b : c);
}
inline void schShrink_w(CvMat *sch, int r, int* pos) {
if (r <= 0) {
printf("Error: r < 0"); throw new cv::Exception();
}
for(int i=0; iheight; ++i)
for(int j=pos[i/r]*r; jwidth - r; ++j)
(CV_MAT_ELEM(*sch, float, i, j)) = CV_MAT_ELEM(*sch, float, i, (j + r));
sch->width = --sch->cols;
}
inline void pixelMove(IplImage *img, CvRect src, CvPoint dst) {
for(int i=src.y, di=dst.y; iheight; ++i, ++di)
for(int j=src.x, dj=dst.x; jwidth; ++j, ++dj)
for(int k=0; knChannels; ++k)
(CV_IMAGE_ELEM(img, uchar, di, dj * img->nChannels + k)) = CV_IMAGE_ELEM(img, uchar, i, j * img->nChannels + k);
}
void seamcarve_wp(IplImage *&img, int h, int reduce_times) {
if (reduce_times <= 0) return;
try {
CvMat* sc = sce;
switch(h) {
case 2:
sc = sch2; break;
case 4:
sc = sch4; break;
case 8:
sc = sch8; break;
case 16:
sc = sch16; break;
default:
break;
}
int *pos = new int[sc->height];
int x, y;
CvMat *sctmp = cvCreateMat(sc->height, sc->width, CV_32F);
for(x = 0; x < sctmp->width; ++x) {
ssetReal(sctmp, 0, x, sgetReal(sc, 0, x));
}
while (reduce_times--) {
for(y = 1; y < sctmp->height; ++y)
for(x = 0; x < sctmp->width; ++x) {
double min = sgetReal(sctmp, y - 1, x);
double n1 = sgetReal(sctmp, y - 1, x - 1, 1e100);
double n2 = sgetReal(sctmp, y - 1, x + 1, 1e100);
min = min3(min, n1, n2);
ssetReal(sctmp, y, x, min + sgetReal(sc, y, x));
}
double m = 1e100; int x0 = 0;
for(x = 0; x < sctmp->width; ++x) {
double n = sgetReal(sctmp, sctmp->height - 1, x);
if (m > n) { m = n; x0 = x ; }
}
// printf("%.4fn", m);
int x = x0;
for(y = sctmp->height - 1; y > 0; --y) {
pos[y] = x;
double val = sgetReal(sctmp, y, x);
double n0 = sgetReal(sctmp, y - 1, x);
double n1 = sgetReal(sctmp, y - 1, x - 1, 1e100);
double n2 = sgetReal(sctmp, y - 1, x + 1, 1e100);
double min = min3(n0, n1, n2);
if (min == n0) ; else if (min == n1) x--; else if (min == n2) x++;
}
pos[0] = x;
// shrink
for(int i=0; iheight; ++i) {
// printf("%d ", pos[i]);
pixelMove(img, cvRect(pos[i]*h + h, i*h, img->width - pos[i] * h, h), cvPoint(pos[i]*h, i*h));
}
//printf("n");
schShrink_w(sctmp, 1, pos);
schShrink_w(sce, h, pos);
if (h >= 2) schShrink_w(sch2, h / 2, pos);
if (h >= 4) schShrink_w(sch4, h / 4, pos);
if (h >= 8) schShrink_w(sch8, h / 8, pos);
if (h >= 16) schShrink_w(sch16, h / 16, pos);
}
cvReleaseMat(≻tmp);
delete[] pos;
} catch (cv::Exception ex) {
printf("%s", ex.what());
}
}
void seamcarve_w(IplImage *&img, int reduce_w) {
// seamcarve_wp(img, 16, reduce_w / 16); reduce_w %= 16;
seamcarve_wp(img, 8, reduce_w / 8); reduce_w %= 8;
seamcarve_wp(img, 4, reduce_w / 4); reduce_w %= 4;
seamcarve_wp(img, 2, reduce_w / 2); reduce_w %= 2;
seamcarve_wp(img, 1, reduce_w / 1);
}
void seamcarve_h(IplImage *&img, int reduce_h) {
}
void seamcarve_calce(IplImage *&src, CvArr *bias, CvMat *sce) {
try {
cvSetZero(sce);
CvMat *sctmp = cvCreateMat(sce->height, sce->width, CV_32F);
cvSobel(src, sctmp, 1, 0); cvSquareAcc(sctmp, sce);
cvSobel(src, sctmp, 0, 1); cvSquareAcc(sctmp, sce);
//cvConvertScale(bias, sctmp, BIAS_WEIGHT);
//cvAdd(sce, sctmp, sce);
cvMul(sce, bias, sce);
// normalization
double min = 0, max = 0;
cvMinMaxLoc(sce, &min, &max, NULL, NULL);
max-=min;
cvConvertScale(sce, sce, 255.0 / max, -min / max * 255.0);
cvReleaseMat(≻tmp);
} catch (cv::Exception ex) {
puts(ex.what());
}
}
void finalize() {
cvReleaseImage(&f);
cvReleaseImage(∈v);
cvReleaseImage(&re);
cvReleaseImage(&im);
cvReleaseImage(&re2);
cvReleaseImage(&im2);
cvReleaseImage(&modular);
cvReleaseImage(&sine);
cvReleaseImage(&cosine);
cvReleaseImage(&logamp);
cvReleaseImage(&resi);
cvReleaseImage(&tmp);
cvReleaseImage(&tmp2);
cvReleaseMat(⪰);
cvReleaseMat(≻h2);
cvReleaseMat(≻h4);
cvReleaseMat(≻h8);
cvReleaseMat(≻h16);
}
CvCapture *capture = 0;
int main(int argc, char* argv[])
{
bool camera = false; int reduce_w = 0, reduce_h = 0;
if (argc != 2) camera = true;
if (argc == 3) {
camera = true;
reduce_w = atoi(argv[1]); reduce_h = atoi(argv[2]);
}
if (argc == 4) {
camera = false;
reduce_w = atoi(argv[2]); reduce_h = atoi(argv[3]);
}
cvNamedWindow(WINDOWNAME,1);
IplImage *img = NULL, *ampl = NULL, *ampls = NULL, *gray = NULL, *grays = NULL;
if (camera) {
capture = cvCaptureFromCAM( 0 );
if (capture == NULL) {
printf("Error: capture failed.n");
return 0;
}
img = cvQueryFrame(capture);
} else {
img = cvLoadImage(argv[1]);
}
gray = cvCreateImage(cvGetSize(img), img->depth, 1);
grays = cvCreateImage(cvSize(SALIENCY_SIZE_W, SALIENCY_SIZE_H), img->depth, 1);
ampl = cvCreateImage(cvGetSize(img), IPL_DEPTH_32F, 1);
ampls = cvCreateImage(cvGetSize(grays), IPL_DEPTH_32F, 1);
init(grays, img);
try {
do {
if (camera)
img = cvQueryFrame(capture);
cvCvtColor(img, gray, CV_RGB2GRAY);
cvResize(gray, grays);
spect_ampl(grays, ampls);
cvResize(ampls, ampl);
// visualization
if (reduce_w > 0 || reduce_h > 0) {
seamcarve_calce(gray, ampl, sce);
cvResize(sce, sch2);
cvResize(sch2, sch4);
cvResize(sch4, sch8);
cvResize(sch8, sch16);
if (reduce_w > 0) seamcarve_w(img, reduce_w);
if (reduce_h > 0) seamcarve_h(img, reduce_h);
cvDrawRect(img, cvPoint(0, 0), cvPoint(img->width - reduce_w, img->height - reduce_h), cvScalar(255));
sce->cols = sce->width = img->width; sce->rows = sce->height = img->height;
sch2->cols = sch2->width = sce->width / 2; sch2->rows = sch2->height = sce->height / 2;
sch4 ->cols = sch4 ->width = sce->width / 4 ; sch4 ->rows = sch4 ->height = sce->height / 4 ;
sch8 ->cols = sch8 ->width = sce->width / 8 ; sch8 ->rows = sch8 ->height = sce->height / 8 ;
sch16->cols = sch16->width = sce->width / 16; sch16->rows = sch16->height = sce->height / 16;
cvShowImage("spec saliency / seam carving energy map", sce);
} else {
int x = 0, y = 0;
float *amplyptr = (float*)(ampl->imageData);
for(char *yptr = img->imageData; y < img->height; ++y, yptr += img->widthStep, amplyptr += ampl->widthStep / sizeof(float)) {
x = 0; float *amplxptr = amplyptr;
for(char *xptr = yptr; x < img->width; ++x, xptr += 3, amplxptr++) {
if (*amplxptr < threshold) {
xptr[0] = xptr[1] = xptr[2] = 255;
}
}
}
cvShowImage("spec saliency / seam carving energy map", ampl);
}
//cvSetImageROI(img, cvRect(0, 0, img->width - reduce_w, img->height - reduce_h));
cvShowImage(WINDOWNAME, img);
if (cvWaitKey(10) == 27) break;
cvResetImageROI(img);
} while(camera);
} catch(cv::Exception ex) {
printf("%sn", ex.what());
}
if (!camera) cvWaitKey();
finalize();
cvDestroyAllWindows();
cvReleaseImage(&img);
cvReleaseImage(&gray);
return 0;
}
