高斯分布的点落入心形曲线的一个解决方案
2017-03-05 13:29
169 查看
给定心形曲线(x2+y2−1)3=x2y3,给定任意一点的坐标(X,Y)其中X~N(X,σx),<
10174
span style="display: inline-block; position: relative; width: 5.369em; height: 0px; font-size: 123%;">Y~N(Y,σy)求点(X,Y)落入心形曲线内的概率。
思路:
以(X,Y)为中心,画出3∗σ半径的椭圆,求和心形曲线相交的体积。注意:心形曲线方程可化为x2+y2−1=x2/3y,满足x2+y2−1<=(x2)1/3y在曲线内。利用心形曲线上下左右都有最大值且约等于正负1。可以设定一个分辨率画出图形。
上代码:
效果图:
心形
p=1.0
p=0.9
p=0.5
p=0.4
10174
span style="display: inline-block; position: relative; width: 5.369em; height: 0px; font-size: 123%;">Y~N(Y,σy)求点(X,Y)落入心形曲线内的概率。
思路:
以(X,Y)为中心,画出3∗σ半径的椭圆,求和心形曲线相交的体积。注意:心形曲线方程可化为x2+y2−1=x2/3y,满足x2+y2−1<=(x2)1/3y在曲线内。利用心形曲线上下左右都有最大值且约等于正负1。可以设定一个分辨率画出图形。
上代码:
# -*- coding: utf-8 -*- import numpy as np import matplotlib.pyplot as plt res=0.01#单位每像素 RES=1/res#像素每单位 block=256 map1=np.ones([block,block]) CX=block/2.0 CY=block/2.0 for y in np.arange(0,block): for x in np.arange(0,block): if (res*(x-CX))**2+(res*(y-CX))**2-1<=(res*np.abs(x-CX))**(2.0/3.0)*(res*(y-CY)): map1[y,x]=0 plt.figure(1) plt.imshow(map1,cmap='gray') sigmax=0.3 sigmay=0.1 X=0 Y=0 l=max(CX+(X-3*sigmax)*RES,0) r=min(CX+(X+3*sigmax)*RES,block) t=max(CY+(Y-3*sigmay)*RES,0) b=min(CY+(Y+3*sigmay)*RES,block) print(l,r,t,b) theta=1/((2.0*np.pi)*(sigmax*sigmay)) ssum=0; for y in np.arange(l,r+1): for x in np.arange(t,b+1): if map1[y,x]==0: map1[y,x]=np.exp(-0.5*((((x-CX)*res-X)/sigmax)**2+(((y-CY)*res-Y)/sigmay)**2)) ssum=ssum+theta*map1[y,x] plt.figure(2) plt.imshow(map1,cmap='gray') #print(ssum/(np.sum(np.sum(f)))) #print(res**2*theta*np.sum(np.sum(f))) print(res**2*ssum) print(np.round(res**2*ssum*10)/10) plt.show()
效果图:
心形
p=1.0
p=0.9
p=0.5
p=0.4
相关文章推荐
- [question] 怎么画一个变量的数据的概率密度分布曲线?
- 人生就是一个高斯分布或者是一个冲击函数
- [原创]《程序员,在路上……》第1节——用OPENEL画出麦克斯维速率分布曲线
- 优化了的过关键点的光滑曲线拟合算法的修正(一个链表的定义)
- 关于在一个form表单里同时上传多个文件和文本信息的解决方案。。。
- 一个简单的解决方案文档-仅做参考.
- 确保只有一个程序实例运行(C#)之解决方案
- 用程序画出麦克斯维速率分布曲线
- Oracle新手最常碰到的6个错误及解决方案-给一个兄弟收集的
- 在XP+vs2002下打开解决方案时提示Web访问失败的一个解决方法
- 一个统计当前在线用户的解决方案
- 一个带复合表头与跨列表项的数据表的DataGrid解决方案
- 莫名其妙的PostBack错误和一个似是而非的解决方案
- 关于上一个接口实现的解决方案
- [总结]C#判断一个string是否可以为数字,五种解决方案!
- iBatis:又一个O/R Mapping解决方案
- 使用displaytag时出现的一个异常及其解决方案
- 一个解决方案对dll的引用要注意的问题
- 如何建立一个灵活的、可配置的导出数据到Excel的解决方案。