用C语言实现瑞利分布，莱斯分布，高斯分布的分布函数

下面共有两个程序，程序2 加入了图形显示
程序1
这个程序就是你要的。
# include "stdioh"
# include "mathh"
# include "stdlibh"
# include "mathh"
# include "dosh"
# define MAX_N 3000 /这个值为N可以定义的最大长度/
# define N 100 /产生随机序列的点数，注意不要大于MAX_N/
/产生均匀分布的随机变量/
void randa(float x,int num);
/产生瑞利分布的随机变量/
void randr(float x,int num);
/产生标准高斯分布的随机变量/
void randn(float x,int num);
/产生莱斯分布的随机变量/
void randl(float x, float a, float b, int num);
void fshow(char name,float x,int num);
main()
{
float x[N];
int i;
/
randa(&x,N);
randr(&x,N);
randl(&x,10,10,N);
/
randn(&x,N);
/此时x[N]就是所需要的高斯分布的序列/
/显示该序列/
fshow("x",&x,N);
getch();
}
void randa(float x,int num)
{
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x[i]=rand();
x[i]=x[i]/32768;
}
}
void randr(float x,int num)
{
float x1[MAX_N];
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x1[i]=rand();
x[i]=x1[i]/32768;
x[i]=sqrt(-2log(x[i]));
}
}
void randn(float x,int num)
{
float x1[MAX_N],x2[MAX_N];
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x1[i]=rand();
x2[i]=rand();
x1[i]=x1[i]/32768;
x2[i]=x2[i]/32768;
x[i]=sqrt(-2log(x1[i]))cos(x2[i]M_PI);
}
}
void randl(float x, float a, float b, int num)
{
float x1[MAX_N],x2[MAX_N];
float temp[MAX_N];
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x1[i]=rand();
x2[i]=rand();
x1[i]=x1[i]/32768;
x2[i]=x2[i]/32768;
temp[i]=sqrt(-2log(x1[i]))cos(x2[i]M_PI);
x2[i]=sqrt(-2log(x1[i]))sin(x2[i]M_PI);
x1[i]=temp[i];
x[i]=sqrt((a+x1[i])(a+x1[i])+(b+x2[i])(b+x2[i]));
}
}
void fshow(char name,float x,int num)
{
int i,sign,L;
float temp;
printf("\n");
printf(name);
printf(" = ");
L=6;
/按照每行6个数据的格式显示/
for(i=0;i<num;i++)
{
temp=i/L;
sign=temp;
if((i-signL)==0) printf("\n");
if(x[i]>0) printf(" %f ",x[i]);
else printf("%f ",x[i]);
}
}
程序 2
以下程序加入了图形显示的效果，因此更加直观，你可以参考一下。
/ 作者 Leo_nanjing
时间 2008510
功能 生成各种分布的随机变量,并显示
/
# include "stdioh"
# include "mathh"
# include "graphicsh"
# include "mathh"
# include "dosh"
# define MAX_N 3000
# define N 1000
void randa(float x,int num);
void randr(float x,int num);
void randn(float x,int num);
void randl(float x, float a, float b, int num);
void fshow(char name,float x,int num);
/用于图形显示的部分/
void init_graphic(unsigned color);
void plotxy(float x, float y, int num,int mode);
void plot(float y,int num, int mode);
float max(float x, int num);
float min(float x, int num);
/画出该随机序列的分布函数曲线/
void plotpdf(float x,int num,int part,int mode);
main()
{
float x[N];
int i;
randn(&x,N);
fshow("x",&x,N);
getch();
/以下为图形显示部分/
init_graphic(0);
/显示随机序列/
plot(&x,N,1);
getch();
/显示其分布函数/
plotpdf(&x,N,20,0);
getch();
}
void randa(float x,int num)
{
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x[i]=rand();
x[i]=x[i]/32768;
}
}
void randr(float x,int num)
{
float x1[MAX_N];
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x1[i]=rand();
x[i]=x1[i]/32768;
x[i]=sqrt(-2log(x[i]));
}
}
void randn(float x,int num)
{
float x1[MAX_N],x2[MAX_N];
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x1[i]=rand();
x2[i]=rand();
x1[i]=x1[i]/32768;
x2[i]=x2[i]/32768;
x[i]=sqrt(-2log(x1[i]))cos(x2[i]M_PI);
}
}
void randl(float x, float a, float b, int num)
{
float x1[MAX_N],x2[MAX_N];
float temp[MAX_N];
int i;
struct time stime;
unsigned seed;
gettime(&stime);
seed=stimeti_hundstimeti_minstimeti_hour;
srand(seed);
for(i=0;i<num;i++)
{
x1[i]=rand();
x2[i]=rand();
x1[i]=x1[i]/32768;
x2[i]=x2[i]/32768;
temp[i]=sqrt(-2log(x1[i]))cos(x2[i]M_PI);
x2[i]=sqrt(-2log(x1[i]))sin(x2[i]M_PI);
x1[i]=temp[i];
x[i]=sqrt((a+x1[i])(a+x1[i])+(b+x2[i])(b+x2[i]));
}
}
void fshow(char name,float x,int num)
{
int i,sign,L;
float temp;
printf("\n");
printf(name);
printf(" = ");
L=6;
for(i=0;i<num;i++)
{
temp=i/L;
sign=temp;
if((i-signL)==0) printf("\n");
if(x[i]>0) printf(" %f ",x[i]);
else printf("%f ",x[i]);
}
}
/以下为图形显示的函数/
void init_graphic(unsigned color)
{
int graphicdriver,graphicmode;
graphicdriver=DETECT;
graphicmode=1;
initgraph(&graphicdriver,&graphicmode,"E:\\turboc2\\");
setbkcolor(color);
}
void plotxy(float x, floaty, int num,int mode)
{
int i;
float max_x,max_y,min_x,min_y;
float x0,y0,x1,y1;
clrscr(0);
cleardevice();
setbkcolor(0);
max_x=max(x,num);
max_y=max(y,num);
min_x=min(x,num);
min_y=min(y,num);
setlinestyle(0,2,1);
line(65,35,65,445);
line(65,445,575,445);
setlinestyle(3,0,1);
line(65,35,575,35);
line(575,35,575,445);
setlinestyle(0,2,1);
if(max_x==min_x)
x0=320;
else
x0=(x[0]-min_x)500/(max_x-min_x)+70;
if(max_y==min_y)
y0=240;
else
y0=480-((y[0]-min_y)400/(max_y-min_y)+40);
if(mode==0) circle(x0,y0,2);
for(i=1;i<num;i++)
{
if(max_x==min_x)
x1=320;
else
x1=(x[i]-min_x)500/(max_x-min_x)+70;
if(max_y==min_y)
y1=240;
else
y1=480-((y[i]-min_y)400/(max_y-min_y)+40);
if(mode==0) circle(x1,y1,2);
line(x0,y0,x1,y1);
x0=x1;y0=y1;
}
printf("\n\n");
printf("%f",max_y);
printf("\n\n\n\n\n\n\n\n\n\n");
printf("\n\n\n");
printf("%f",(max_y+min_y)/2);
printf("\n\n\n\n\n\n\n\n\n\n");
printf("\n\n");
printf("%f",min_y);
printf("\n %f",min_x);
printf(" ");
printf("%f",(max_x+min_x)/2);
printf(" ");
printf("%f",max_x);
}
void plot(floaty, int num,int mode)
{
int i;
float max_x,max_y,min_x,min_y;
float x0,y0,x1,y1;
float x[MAX_N];
clrscr(0);
cleardevice();
setbkcolor(0);
for(i=0;i<num;i++) x[i]=i+1;
max_x=max(x,num);
max_y=max(y,num);
min_x=min(x,num);
min_y=min(y,num);
setlinestyle(0,2,1);
line(65,35,65,445);
line(65,445,575,445);
setlinestyle(3,0,1);
line(65,35,575,35);
line(575,35,575,445);
setlinestyle(0,2,1);
if(max_x==min_x)
x0=320;
else
x0=(x[0]-min_x)500/(max_x-min_x)+70;
if(max_y==min_y)
y0=240;
else
y0=480-((y[0]-min_y)400/(max_y-min_y)+40);
if(mode==0) circle(x0,y0,2);
for(i=1;i<num;i++)
{
if(max_x==min_x)
x1=320;
else
x1=(x[i]-min_x)500/(max_x-min_x)+70;
if(max_y==min_y)
y1=240;
else
y1=480-((y[i]-min_y)400/(max_y-min_y)+40);
if(mode==0) circle(x1,y1,2);
line(x0,y0,x1,y1);
x0=x1;y0=y1;
}
printf("\n\n");
printf("%f",max_y);
printf("\n\n\n\n\n\n\n\n\n\n");
printf("\n\n\n");
printf("%f",(max_y+min_y)/2);
printf("\n\n\n\n\n\n\n\n\n\n");
printf("\n\n");
printf("%f",min_y);
printf("\n %f",min_x);
printf(" ");
printf("%f",(max_x+min_x)/2);
printf(" ");
printf("%f",max_x);
}
void plotpdf(float x,int num,int part,int mode)
{
int i,j;
float max_x,min_x,round,deltax,up,down,sum;
float xl[MAX_N],yl[MAX_N];
sum=0;
max_x=max(x,num);
min_x=min(x,num);
round=max_x-min_x;
deltax=round/part;
xl[0]=min_x;
for(i=1;i<=part;i++)
{
xl[i]=min_x+deltaxi;
yl[i-1]=0;
up=xl[i];
down=xl[i-1];
for(j=0;j<num;j++)
{
if((x[j]<up) && (x[j]>=down)) yl[i-1]=yl[i-1]+1;
}
yl[i-1]=yl[i-1]/num/deltax;
}
for(i=0;i<part;i++) sum=sum+yl[i];
plotxy(&xl,&yl,part,mode);
}
float max(float x, int num)
{
int i;
float max;
max=x[0];
for(i=1;i<num;i++)
{
if(x[i]>max) max=x[i];
}
return max;
}
float min(float x, int num)
{
int i;
float min;
min=x[0];
for(i=1;i<num;i++)
{
if(x[i]<min) min=x[i];
}
return min;
}

一般加入白噪声，只得是加入高斯白噪声 利用randn产生高斯分布噪声加到原来的图像上就可以了 假设你原来的图像的数据存在矩阵变量a中 那么 noise=Arandn(size(a)); 其中A是需要加噪声的方差值，A的取值确定噪声的幅度大小 一般比图像信号幅度

---