博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
计算几何知识
阅读量:3951 次
发布时间:2019-05-24

本文共 6150 字,大约阅读时间需要 20 分钟。

文章目录

注意

  • 浮点数0,有可能输出-0。需要x = sgn(x) == 0 ? 0 : x;处理一下。
  • 精度要尽量开高一点,一般1e-10比较合适。
  • double atan2(double y,double x)反正切函数,能判断象限输出 ( − π , π ] (-\pi,\pi] (π,π]范围的角度,且避免除零问题。

向量叉乘的意义

方向:右手定则

大小:两向量组成的平行四边形的面积

直线和线段是否相交

跨立实验,judge_point其中一个为0,或者异号。

struct Point{
double x, y; double operator* (const Point &p) const {
return x*p.y - p.x*y; } Point operator- (const Point &p) const {
return Point{
x-p.x, y-p.y}; } double len(){
return sqrt(x*x+y*y); }};struct Line{
Point a,b; Point dir() const{
return a - b; } double judge_point(const Point &p) const{
return this->dir() * (a - p); }};

两线段是否相

在这里插入图片描述

bool iscross(Line s1,Line s2){
return min(s1.a.x,s1.b.x) <= max(s2.a.x,s2.b.x) && min(s1.a.y,s1.b.y) <= max(s2.a.y,s2.b.y) && min(s2.a.x,s2.b.x) <= max(s1.a.x,s1.b.x) && min(s2.a.y,s2.b.y) <= max(s1.a.y,s1.b.y) && s1.judge_point(s2.a) * s1.judge_point(s2.b) < 1e-8 && s2.judge_point(s1.a) * s2.judge_point(s1.b) < 1e-8;}

求两直线交点

先判断两直线是否平行(方向向量叉积为0)或共线(并有一个点)。再用以下公式求

x0 = ((x3-x4) * (x2*y1 - x1*y2) - (x1-x2) * (x4*y3 - x3*y4)) / ((x3-x4) * (y1-y2) - (x1-x2) * (y3-y4));y0 = ((y3-y4) * (y2*x1 - y1*x2) - (y1-y2) * (y4*x3 - y3*x4)) / ((y3-y4) * (x1-x2) - (y1-y2) * (x3-x4));

坐标系旋转

在平面坐标上,任意点P(x1,y1),绕一个坐标点Q(x2,y2)旋转 θ θ θ角度后,新的坐标设为(x, y)的计算公式:

x= (x1 - x2)*cos(θ) - (y1 - y2)*sin(θ) + x2 ;y= (x1 - x2)*sin(θ) + (y1 - y2)*cos(θ) + y2 ;

(原理:先把坐标系中心平移到Q点, 把P点在极坐标下旋转 θ θ θ, 再把坐标轴平移回去)

极角排序

通过叉乘和距离判断

bool cmp(const Point &a, const Point &b){
Point da = a - mark; Point db = b - mark; switch(sgn(da ^ db)){
case 1: return 1; case -1: return 0; default: return da.len() < db.len(); }}

判断一个点是否在多边形内

回转数法

double get_angle(Point a, Point b, Point c){
a = a - c; b = b - c; return acos((a^b) / (a.len() * b.len())); // ^是内积}bool is_in_ploygn(Point p){
// 回转数法 提前把v[0]放到最后面 double tot = 0; for(int i=0;i
0) tot += get_angle(v[i],v[i+1],p); else tot -= get_angle(v[i],v[i+1],p); } tot = fabs(tot); if(sgn(tot) == 0)return false; // out if(sgn(tot - 2 * PI) == 0) return true; // in if(sgn(tot - PI) == 0) return true; // 在多边形边上 不在拐角点上 else return true; // 在多边形拐点上}

判断一个多边形是否为凸(已考虑顺逆时针)

bool is_convex(){
// 提前先把前两个点加到v后面 int f=0; for(int tmp, i=0;i

Pick定理

以整点坐标组成的多边形面积为

在这里插入图片描述
内部格点数目 i i i、边上格点数目 b b b

多边形面积

任意一个多边形的面积等于按顺序求相邻两个点与原点组成的向量的叉积之和除2的绝对值(注意正负)。

二维凸包

Graham 扫描法

  1. 找出最左下角的点,作为凸包第一个点
  2. 以这个点为key,对剩下的点进行极角排序(画图可知,凸包中含的点序列必然在此时升序,且第一个点和最后一个点必然为凸包上点)
  3. 依次枚举,判断当前点和凸包前两个点是否为右拐,右拐则前一个点不是凸包上的点,退栈。直到可以左拐或栈只有一个点为止。加入该点。
    凸包示意图
int id = 0;for(int i=1;i
0 && sgn( (sta[top] - sta[top-1]) * (v[i] - sta[top]) ) < 0) top--; // 踢掉右转的 sta[++top] = v[i];}

上凸包,下凸包

有时候,需要单求上凸包,或下凸包。求凸包方法,与上面类似,但是一开始不是用极角排序,而是按第一关键字x,第二关键字y排序。

bool cmp(const Point &a, const Point &b){
return a.x < b.x || a.x == b.x && a.y < b.y;}

下凸包,即保证叉乘大于零

sta[top = 1] = p[1];for(int i = 2; i <= ln; ++i){
// 下凸包 while(top > 1 && ((sta[top] - sta[top-1]) ^ (p[i] - sta[top])) <= 0)--top; sta[++top] = p[i];}

上凸包,即保证叉乘小于零

sta[top = 1] = p[1];for(int i = 2; i <= ln; ++i){
// 上凸包 while(top > 1 && ((sta[top] - sta[top-1]) ^ (p[i] - sta[top])) >= 0)--top; sta[++top] = p[i];}

旋转卡壳

在这里插入图片描述

在这里插入图片描述

被两条线夹逼着整个凸包,就产生了一对对踵点。可以证明对踵点的个数不超过 3 N / 2 3N/2 3N/2个 也就是说对踵点的个数 O ( N ) O(N) O(N)的。旋转卡壳就是用来求对踵点对的 O ( n ) O(n) O(n)算法。由于,固定一条边后,求离这条边最远的点,是一个单峰函数。并且容易看出,当边逆时针转动时,最远的点必然也在逆时针转动。就可得到以下代码,通过叉积得到三角形的面积,间接比较边到点之间的距离。

int area(int i, int j, int k){
int res = (sta[j] - sta[i]) ^ (sta[k] - sta[i]); return res;}
key = a[0]; //前面已保证a[0]是y最小,x最小的点    sort(a + 1, a + n, cmp); // 极角排序    int top = 0;    sta[0] = a[0]; // 求凸包    for(int i = 1; i < n; ++i){
while(top > 0 && ((sta[top] - sta[top - 1]) ^ (a[i] - sta[top])) <= 0) top--; sta[++top] = a[i]; } int ans = 0; sta[++top] = sta[0]; int j = 1; // 旋转卡壳 for(int i = 0; i < top; ++i){
while(area(i, i + 1, j) < area(i, i + 1, j + 1) ){
j = (j + 1) % top; } ans = max(ans, (sta[i] - sta[j]).len()); }

扫描线

在这里插入图片描述

求矩形组成的面积或周长。用线段树统计目前每一段的计数,计数大于0则统计该段长度。

#include 
#define endl '\n'#define IOS ios::sync_with_stdio(false); cin.tie(nullptr), cout.tie(nullptr);using namespace std;const int N = 2e2 + 5;double valx[N];struct Point{
double x1, x2, y; int flag; bool operator< (const Point &ot) const &{
return y < ot.y; }}p[N];#define ls rt << 1#define rs rt << 1 | 1double sum[N << 3];int lazy[N << 3];void pushup(int rt, int l, int r){
if(lazy[rt]){
sum[rt] = valx[r - 1] - valx[l - 1]; }else if(r + 1 != l)sum[rt] = sum[ls] + sum[rs]; else sum[rt] = 0;}void update(int rt, int l, int r, int lx, int rx, int f){
if(lx <= l && r <= rx){
lazy[rt] += f; pushup(rt, l, r); return; } int mid = l + r >> 1; if(lx < mid)update(ls, l, mid, lx, rx, f); if(rx > mid)update(rs, mid, r, lx, rx, f); pushup(rt, l, r);}int getpos(double x, int n){
return lower_bound(valx, valx + n, x) - valx + 1;}int main(){
IOS; int n, cas = 1; while(cin >> n && n){
for(int i = 0; i < n; ++i){
double x1, x2, y1, y2; cin >> x1 >> y1 >> x2 >> y2; valx[i] = x1, valx[i + n] = x2; p[i].x1 = x1, p[i].x2 = x2, p[i].y = y1, p[i].flag = 1; p[n + i].x1 = x1, p[n + i].x2 = x2, p[n + i].y = y2, p[n + i].flag = -1; } memset(sum, 0, sizeof(sum)); memset(lazy, 0, sizeof(lazy)); sort(valx, valx + 2 * n); int m = unique(valx, valx + 2 * n) - valx; sort(p, p + 2 * n); double ans = 0; for(int i = 0; i < 2 * n; ++i){
//cout << p[i].y << " " << p[i].x1 << " " << p[i].x2 << " " << p[i].flag << " sum=" << sum[1] << endl; if(i)ans += (p[i].y - p[i - 1].y) * sum[1]; update(1, 1, m, getpos(p[i].x1, m), getpos(p[i].x2, m), p[i].flag); } cout << "Test case #" << cas++ << "\nTotal explored area: " << fixed << setprecision(2) << ans << endl << endl; }}

转载地址:http://zwuzi.baihongyu.com/

你可能感兴趣的文章
Spring学习之Filter、Interceptor、Aop实现与区别
查看>>
Spring 添加@Autowired注释, 注入对象却为空
查看>>
springSecurity学习
查看>>
通过Java的api操作redis
查看>>
jquery基本选择器
查看>>
linux学习之shell字符串大小写转换
查看>>
Linux下用base64对字符串进行加密解密
查看>>
H5走迷宫小游戏
查看>>
mysql建表 表名与关键字冲突
查看>>
mysql 创建单表外键关联多表
查看>>
postman使用
查看>>
ClassNotFoundException和NoClassDefFoundError的区别
查看>>
Tomcat Connector三种运行模式(BIO, NIO, APR)的比较和优化
查看>>
spring注解@Primary与@Qualifier
查看>>
annotation之@Autowired、@Inject、@Resource三者区别
查看>>
idea启动微服务找不到配置文件
查看>>
Java通过反射机制调用某个类的方法
查看>>
字节跳到面试题
查看>>
Linux查看物理CPU个数
查看>>
Linux学习之网络IO,磁盘io
查看>>