Z-Ordering 编码及查询c++实现(简称Z2)

网上搜了很多Z-Ordering实现没搜到,通过 sfcurve-master和geomesa-geomesa-3.2.2 得scala代码改编而来, 环境为C++, vs2015, 理论上windows和Linux都可以用.

#include "stdafx.h"
#include <vector>
#include <list>
#include <queue>
#include <algorithm>
#include <iostream>
#include <deque>
#include <math.h>
#include <time.h>
#include <functional>  
#include <memory>
#include <set>
#include <map>
#include <queue>
#include <string>
#include <memory>
#include <functional>
#include <math.h>
using namespace std;



typedef unsigned long long ull;

struct ZRange {
	ull min = 0, max = 0, mid = 0, length = 0;

	ZRange(ull umin, ull umax)
	{
		min = umin;
		max = umax;
		mid = (max + min) >> 1;
		length = max - min + 1;
	}
	bool contains(ull bits)
	{
		return   bits >= min && bits <= max;
	}
	bool contains(const ZRange& r)
	{
		return contains(r.min) && contains(r.max);
	}
	bool overlaps(const ZRange& r) {
		return  contains(r.min) || contains(r.max);
	}
};


struct Z2
{
	int Dimensions = 2;
	int BitsPerDimension = 31;
	int TotalBits = 62;
	ull MaxMask = 0x7fffffffL;
	ull z = 0;
	int d0 = 0, d1 = 0;
	int DefaultRecurse = 7;
	Z2(ull Z) :z(Z)
	{
		d0 = dim(0);
		d1 = dim(1);
	}
	Z2(int x, int y)
	{
		z = split(x) | (split(y) << 1);
		d0 = dim(0);
		d1 = dim(1);
	}
	//def <(other: Z2) = z < other.z
	//def >(other: Z2) = z > other.z
	//def <= (other: Z2) = z <= other.z
	//def >= (other: Z2) = z >= other.z
	//def + (offset: Long) = new Z2(z + offset)
	//def - (offset: Long) = new Z2(z - offset)
	//def == (other: Z2) = other.z == z

	static int combine(ull z)
	{
		ull x = z & 0x5555555555555555L;
		x = (x ^ (x >> 1)) & 0x3333333333333333L;
		x = (x ^ (x >> 2)) & 0x0f0f0f0f0f0f0f0fL;
		x = (x ^ (x >> 4)) & 0x00ff00ff00ff00ffL;
		x = (x ^ (x >> 8)) & 0x0000ffff0000ffffL;
		x = (x ^ (x >> 16)) & 0x00000000ffffffffL;
		return x;
	}
	ull split(ull value)
	{
		ull x = value & MaxMask;
		x = (x ^ (x << 32)) & 0x00000000ffffffffL;
		x = (x ^ (x << 16)) & 0x0000ffff0000ffffL;
		x = (x ^ (x << 8)) & 0x00ff00ff00ff00ffL; // 11111111000000001111111100000000..
		x = (x ^ (x << 4)) & 0x0f0f0f0f0f0f0f0fL; // 1111000011110000
		x = (x ^ (x << 2)) & 0x3333333333333333L; // 11001100..
		x = (x ^ (x << 1)) & 0x5555555555555555L;// 1010...
		return x;
	}
	void decode(int &x, int &y)
	{
		x = combine(z);
		y = combine(z >> 1);
	}
	int dim(int i) {
		return   Z2::combine(z >> i);
	}
	Z2	  mid(Z2& p)
	{
		int x, y;
		decode(x, y);
		int px, py;
		p.decode(px, py);
		return   Z2((x + px) >> 1, (y + py) >> 1); // overflow safe mean
	}
	//包含
	static bool  contains(ZRange & range, ull value)
	{
		Z2 z2tmp(value);
		int x, y;
		z2tmp.decode(x, y);
		return x >= Z2(range.min).d0 && x <= Z2(range.max).d0 && y >= Z2(range.min).d1 && y <= Z2(range.max).d1;
	}
	static bool overlaps(int a1, int a2, int b1, int b2)
	{
		return  max(a1, b1) <= min(a2, b2);
	}
	static bool  overlaps(ZRange& range, ZRange& value)
	{
		return	overlaps(Z2(range.min).d0, Z2(range.max).d0, Z2(value.min).d0, Z2(value.max).d0) &&
			overlaps(Z2(range.min).d1, Z2(range.max).d1, Z2(value.min).d1, Z2(value.max).d1);
	}
};
struct BitNormalizedDimension {
	double min, max;
	int maxIndex;
	double normalizer;
	double denormalizer;
	BitNormalizedDimension()
	{

	}
	BitNormalizedDimension(const BitNormalizedDimension & o)
	{
		min = o.min;
		max = o.max;
		maxIndex = o.maxIndex;
		normalizer = o.denormalizer;
		denormalizer = o.denormalizer;
	}
	BitNormalizedDimension(double dmin, double dmax, int precision)
	{
		if (!(precision > 0 && precision < 32))
			return;
		min = dmin;
		max = dmax;
		ull bins = pow(2, precision);// .toLong1 << precision;
		normalizer = bins / (max - min);
		denormalizer = (max - min) / bins;
		maxIndex = (bins - 1);
	}
	int normalize(double x) {
		if (x >= max)
			return maxIndex;
		else
			return  floor((x - min) * normalizer);
	}
	double denormalize(int x) {
		if (x >= maxIndex)
		{
			return min + (maxIndex + 0.5) * denormalizer;
		}
		else
		{
			return min + (x + 0.5) * denormalizer;
		}
	}
};
struct ZOrderingSpaceFillingCurve
{
	BitNormalizedDimension lon;
	BitNormalizedDimension lat;
	ZOrderingSpaceFillingCurve()
	{

	}
	ZOrderingSpaceFillingCurve(int precision)
	{
		lon = BitNormalizedDimension(-180., 180., precision);
		lat = BitNormalizedDimension(-180.0, 180.0, precision);
	}
	//正算
	ull  index(double x, double y, bool lenient = false)
	{
		if (lenient)
		{
			return lenientIndex(x, y);
		}
		return Z2(lon.normalize(x), lat.normalize(y)).z;
	}
	//正算
	ull lenientIndex(double x, double y)
	{
		double bx = x;
		if (x < lon.min) { x = lon.min; }
		else if (x > lon.max) { x = lon.max; }

		double by = y;
		if (y < lat.min) { y = lat.min; }
		else if (y > lat.max) { y = lat.max; }

		return Z2(lon.normalize(bx), lat.normalize(by)).z;
	}
	//反算
	void invert(long long z, double& x, double & y)
	{
		int xr, yc;
		Z2(z).decode(xr, yc);
		x = lon.denormalize(xr);
		y = lat.denormalize(yc);
	}

	List<ZRange> ranges(double xmin, double ymin, double xmax, double ymax,
		int precision,
		int maxRanges)
	{
		//val zbounds = xy.map{ case (xmin, ymin, xmax, ymax) = > ZRange(index(xmin, ymin), index(xmax, ymax)) }
		//	Z2.zranges(zbounds.toArray, precision, maxRanges)
       // 这里没实现, 其实就是遍历格网,然后全部取出来
	}
};

int main()
{
	ZOrderingSpaceFillingCurve f(31);
	std::vector<double> coords = { 78,15,112,45 };
	ull g = f.index(78, 15);
	double x, y;
	//正反算测试
	f.invert(g, x, y);
	g=f.index(90, 19);
	f.invert(g, x, y);
	g=f.index(115, 54);
	f.invert(g, x, y);
	g=f.index(75, 19);
	f.invert(g, x, y);
	g=f.index(79, 13);
	f.invert(g, x, y);
	g=f.index(79, 46);
	f.invert(g, x, y);

	ZRange r1(f.index(78, 15), f.index(112, 45));
	ZRange rin(f.index(78, 15), f.index(112, 45));
	ZRange rover(f.index(90, 19), f.index(115, 54));
	ZRange rout(f.index(113, 15), f.index(118, 56));
	//点测试
	bool tt = Z2::contains(r1, f.index(90, 19));
	bool t1f = Z2::contains(r1, f.index(115, 54));
	bool t2f = Z2::contains(r1, f.index(75, 19));
	bool t3f = Z2::contains(r1, f.index(79, 13));
	bool t4f = Z2::contains(r1, f.index(79, 46));
	//矩形测试 
	bool true1 = Z2::overlaps(r1,rin);
	bool true2 = Z2::overlaps(r1,rover);
	bool false1 = Z2::overlaps(r1,rout);
	return 0;
}