golang求多边形相交面积

直接上代码, 根据网上c++版本翻译

package main

import (
	"fmt"
	"math"
)

const maXn int = 300
const eps float64 = 1e-6

//位置标识
func dcmp(X float64) int {
	if X > eps {
		return 1
	} else if X < -eps {
		return -1
	} else {
		return 0
	}
}

type Point struct {
	X, Y float64
}

func cross(a, b, c Point) float64 {
	return (a.X-c.X)*(b.Y-c.Y) - (b.X-c.X)*(a.Y-c.Y) //叉积公式
}

//计算线段ab和cd的交点坐标
func intersection(a, b, c, d Point) Point {
	var p = a
	t := ((a.X-c.X)*(c.Y-d.Y) - (a.Y-c.Y)*(c.X-d.X)) / ((a.X-b.X)*(c.Y-d.Y) - (a.Y-b.Y)*(c.X-d.X))
	p.X += (b.X - a.X) * t
	p.Y += (b.Y - a.Y) * t
	return p
}

//计算多边形面积,将多边形拆解成连续三个顶点组合成的多个三角形进行计算,这个循环计算一次其实是计算两次多边形的面积。
func PolygonArea(p []Point, n int) float64 {
	if n < 3 {
		return 0.0
	}
	s := p[0].Y * (p[n-1].X - p[1].X)
	for i := 1; i < n-1; i++ {
		s += p[i].Y * (p[i-1].X - p[i+1].X)
	}
	s += p[n-1].Y * (p[n-2].X - p[0].X)
	return math.Abs(s * 0.5)
}

func CPIA(a, b []Point, na, nb int) float64 { //ConveXPolYgonIntersectArea
	p := make([]Point, 20)
	tmp := make([]Point, 20)
	var tn, sflag, eflag int
	//eflag = -1
	copy(p, b[:nb])
	//p = append(p[0:], b[:nb]...)
	for i := 0; i < na && nb > 2; i++ {
		if i == na-1 {
			sflag = dcmp(cross(a[0], p[0], a[i]))
		} else {
			sflag = dcmp(cross(a[i+1], p[0], a[i]))
		}
		tn = 0
		for j := tn; j < nb; j++ {
			if sflag >= 0 {
				tmp[tn] = p[j]
				tn++
			}
			if i == na-1 {
				if j == nb-1 {
					eflag = dcmp(cross(a[0], p[0], a[i]))
				} else {
					eflag = dcmp(cross(a[0], p[j+1], a[i])) //计算下一个连续点在矢量线段的位置
				}
			} else {
				if j == nb-1 {
					eflag = dcmp(cross(a[i+1], p[0], a[i]))
				} else {
					eflag = dcmp(cross(a[i+1], p[j+1], a[i]))
				}
			}
			if (sflag ^ eflag) == -2 { //1和-1的异或为-2,也就是两个点分别在矢量线段的两侧
				if i == na-1 {
					if j == nb-1 {
						tmp[tn] = intersection(a[i], a[0], p[j], p[0]) //求交点
						tn++
					} else {
						tmp[tn] = intersection(a[i], a[0], p[j], p[j+1])
						tn++
					}
				} else {
					if j == nb-1 {
						tmp[tn] = intersection(a[i], a[i+1], p[j], p[0])
						tn++
					} else {
						tmp[tn] = intersection(a[i], a[i+1], p[j], p[j+1])
						tn++
					}
				}
			}
			sflag = eflag
		}
		copy(p, tmp)
		nb = tn
		p[nb] = p[0]
	}
	if nb < 3 {
		return 0.0
	}
	return PolygonArea(p, nb)
}
func swap(a, b *Point) {
	*a, *b = *b, *a
}
func SPIA(a, b []Point, na, nb int) float64 { 
	var i, j int
	t1 := make([]Point, na)
	t2 := make([]Point, nb)
	var res float64 = 0
	var num1, num2 float64
	t1[0] = a[0]
	t2[0] = b[0]
	for i = 2; i < na; i++ {
		t1[1] = a[i-1]
		t1[2] = a[i]
		num1 = float64(dcmp(cross(t1[1], t1[2], t1[0]))) //根据差积公式来计算t1[2]在矢量线段(t1[0], t1[1])的左侧还是右侧,
		//值为负数在矢量线段左侧,值为正数在矢量线段右侧
		if num1 < 0 {
			swap(&t1[1], &t1[2]) // 按逆时针进行排序
		}
		for j = 2; j < nb; j++ {
			t2[1] = b[j-1]
			t2[2] = b[j]
			num2 = float64(dcmp(cross(t2[1], t2[2], t2[0])))
			if num2 < 0 {
				swap(&t2[1], &t2[2])
			}
			res += CPIA(t1, t2, 3, 3) * num1 * num2
		}
	}
	return res
}

func main() {
	var p1, p2 []Point
	var n1, n2 int
	//以四边形为例
	p1 = []Point{{X: 0, Y: 0},
		{0, 5},
		{5, 5},
		{5, 0},
	}
	p2 = []Point{{X: 0, Y: 0},
		{0, 5},
		{5, 5},
		{4, 10},
	}
	n1 = len(p1)
	n2 = len(p2)
	A1 := PolygonArea(p1, n1)
	A2 := PolygonArea(p2, n2)
	crossArea := SPIA(p1, p2, n1, n2)
	fmt.Println("A1 面积:", A1)
	fmt.Println("A2 面积:", A2)
	fmt.Println("相交面积 crossArea:", crossArea)
}

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值