svmsvm.cpp

变化检测源程序

// svmsvm.cpp : Defines the entry point for the console application.
//
#include "svm.h"
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <ctype.h>

#include "stdafx.h"
#define Malloc(type,n) (type *)malloc((n)*sizeof(type))
#define  WIDTH 496
#define  HEIGHT 496
#define BandNum 6
#define BYTE unsigned char
void read_problem2(const char *filename);

void read_problem(const char *filename);

struct svm_parameter param;		// set by parse_command_line
struct svm_problem prob;		// set by read_problem

struct svm_problem prob2;		// set by read_problem
struct svm_model *model;
struct svm_node *x_space;

int nr_fold;
struct svm_node *x;
struct svm_node *x2;

int main(int argc, char **argv)
{
	char input_file_name[]="e:\\windows\\traina.txt";
	char input_file_name2[]="e:\\windows\\trainb.txt";

	char model_file_name[]="e:\\windows\\trainab.model";
	const char strbase[]="e:\\changetest4\\98newbip";
	const char strchange[]="e:\\changetest4\\20newbip";
	//const char strbase[]="e:\\changetest4\\1998subbip496";
	//const char strchange[]="e:\\changetest4\\2000subbip496";

	const char stryanmu[]="e:\\changetest4\\yanmuall.raw";
	const char *error_msg;


	// default values
	param.svm_type = ONE_CLASS_MY;
	param.kernel_type = MY;
	param.degree = 3;
	param.gamma =0.1;
	param.coef0 = 0;
	param.nu =0.001;
	param.cache_size = 40;
	param.C = 1;
	param.eps = 1e-3;
	param.p = 0.1;
	param.shrinking = 1;
	param.probability = 0;
	param.nr_weight = 0;
	param.weight_label = NULL;
	param.weight = NULL;

	read_problem(input_file_name);
	read_problem2(input_file_name2);

	error_msg = svm_check_parameter(&prob,&param);
	if(error_msg)
	{
		fprintf(stderr,"Error: %s\n",error_msg);
		exit(1);
	}
	error_msg = svm_check_parameter(&prob2,&param);
	if(error_msg)
	{
		fprintf(stderr,"Error: %s\n",error_msg);
		exit(1);
	}

	FILE *fp_inbase,*fp_inchange,*fp_yanmu;
	fp_inbase=fopen(strbase,"rb");
	fp_inchange=fopen(strchange,"rb");

	fp_yanmu=fopen(stryanmu,"rb");
	int nsize=WIDTH*HEIGHT;
	unsigned char* pchange=new unsigned char[nsize*BandNum+1];
	unsigned char* pbase=new unsigned char[nsize*BandNum+1];
	unsigned char* pyanmu=new unsigned char[nsize+1];
	if(pbase==NULL||pchange==NULL||pyanmu==NULL)
		printf("内存不足！\n");
	fread(pbase,1,nsize*BandNum,fp_inbase);//读入数据
	fread(pchange,1,nsize*BandNum,fp_inchange);
	fread(pyanmu,1,nsize,fp_yanmu);
	int i,j,k;
	BYTE **bimage=new BYTE*[HEIGHT];
	for(i=0;i<HEIGHT;i++)
		bimage[i]=pbase+i*BandNum*WIDTH;
	BYTE **aimage=new BYTE*[HEIGHT];
	for(i=0;i<HEIGHT;i++)
		aimage[i]=pchange+i*BandNum*WIDTH;

	x=new struct svm_node[7];
	x2=new struct svm_node[7];


	unsigned char*yanmu=new unsigned char[nsize];
	memset(yanmu,255,nsize);

	model = svm_trainmy(&prob,&prob2,&param);
	double predict=0.0;
	int correct=0;
	int total=0;
	for (i=0;i<HEIGHT;i++)
	{
		for (j=0;j<WIDTH;j++)
		{
			if (pyanmu[i*WIDTH+j]==0)
			{
				for(k=0;k<BandNum;k++)
				{
					x[k].index=k+1;
					x[k].value=(double)aimage[i][j*BandNum+k]/256;
					x2[k].index=k+1;
					x2[k].value=(double)bimage[i][j*BandNum+k]/256;
				}
				x[6].index=-1;
				x2[6].index=-1;
				predict=svm_predict(model,x,x2);
				if(predict>0) 
				{
					correct++;
					yanmu[i*WIDTH+j]=0;
				}
				total++;
			}
			
		}
	}

	FILE *pResult;
	pResult=fopen("c:\\temp\\c.raw","wb");
    fwrite(yanmu,sizeof(unsigned char),nsize,pResult);
	// 关闭文件
	fclose(pResult);
	delete []yanmu;

	svm_save_model(model_file_name,model);
	svm_destroy_model(model);
	svm_destroy_param(&param);

	printf("Accuracy = %g%% (%d/%d)\n",
		(double)correct/total*100,correct,total);

	free(x);
	free(x2);
	free(prob.y);
	free(prob.x);
	free(prob2.y);//
	free(prob2.x);//
	free(x_space);
	fclose(fp_inbase);//=fopen(strbase,"rb");
	fclose(fp_inchange);//=fopen(strchange,"rb");
	fclose(fp_yanmu);//=fopen(stryanmu,"rb");
	delete []pchange;
	delete []pbase;
	delete []pyanmu;
	delete []bimage;
	delete []aimage;

	return 0;
}



// read in a problem (in svmlight format)

void read_problem(const char *filename)//读数据
{
	int elements, max_index, i, j;
	FILE *fp = fopen(filename,"r");
	
	if(fp == NULL)
	{
		fprintf(stderr,"can't open input file %s\n",filename);
		exit(1);
	}

	prob.l = 0;//样本个数
	elements = 0;//每个样本的维数
	while(1)
	{
		int c = fgetc(fp);
		switch(c)
		{
			case '\n':
				++prob.l;
				// fall through,
				// count the '-1' element
			case ':':
				++elements;
				break;
			case EOF:
				goto out;
			default:
				;
		}
	}
out:
	rewind(fp);

	prob.y = Malloc(double,prob.l);
	prob.x = Malloc(struct svm_node *,prob.l);
	x_space = Malloc(struct svm_node,elements);

	max_index = 0;
	j=0;
	for(i=0;i<prob.l;i++)
	{
		double label;
		prob.x[i] = &x_space[j];
		fscanf(fp,"%lf",&label);
		prob.y[i] = label;
		while(1)
		{
			int c;
			do {
				c = getc(fp);
				if(c=='\n') goto out2;
			} while(isspace(c));
			ungetc(c,fp);
			fscanf(fp,"%d:%lf",&(x_space[j].index),&(x_space[j].value));
			x_space[j].value/=256.0;
			++j;
		}	
out2:
		if(j>=1 && x_space[j-1].index > max_index)
			max_index = x_space[j-1].index;
		x_space[j++].index = -1;
	}

	if(param.gamma == 0)
		param.gamma = 1.0/max_index;

	fclose(fp);
}

void read_problem2(const char *filename)//读数据
{
	int elements, max_index, i, j;
	FILE *fp = fopen(filename,"r");
	
	if(fp == NULL)
	{
		fprintf(stderr,"can't open input file %s\n",filename);
		exit(1);
	}

	prob2.l = 0;
	elements = 0;
	while(1)
	{
		int c = fgetc(fp);
		switch(c)
		{
			case '\n':
				++prob2.l;
				// fall through,
				// count the '-1' element
			case ':':
				++elements;
				break;
			case EOF:
				goto out;
			default:
				;
		}
	}
out:
	rewind(fp);

	prob2.y = Malloc(double,prob2.l);
	prob2.x = Malloc(struct svm_node *,prob2.l);
	x_space = Malloc(struct svm_node,elements);

	max_index = 0;
	j=0;
	for(i=0;i<prob2.l;i++)
	{
		double label;
		prob2.x[i] = &x_space[j];
		fscanf(fp,"%lf",&label);
		prob2.y[i] = label;
		while(1)
		{
			int c;
			do {
				c = getc(fp);
				if(c=='\n') goto out2;
			} while(isspace(c));
			ungetc(c,fp);
			fscanf(fp,"%d:%lf",&(x_space[j].index),&(x_space[j].value));
			x_space[j].value/=256.0;
			++j;
		}	
out2:
		if(j>=1 && x_space[j-1].index > max_index)
			max_index = x_space[j-1].index;
		x_space[j++].index = -1;
	}

	if(param.gamma == 0)
		param.gamma = 1.0/max_index;

	fclose(fp);
}