通過雷達、激光掃描、立體攝像機等三維測量設備獲取的點雲數據,具有數據量大、分布不均勻等特點。作為三維領域中一個重要的數據來源,點雲數據主要是表征目標表面的海量點集合,並不具備傳統網格數據的集合拓撲信息。所以點雲數據處理中最為核心的問題就是建立離散點間的拓撲關系,實現基於鄰域關系的快速查找。
建立空間索引在點雲數據處理中已被廣泛應用,常見空間索引一般是自頂向下逐級划分空間的各種空間索引結構,比較有代表性的包括 BSP樹、 KD樹、 KDB樹、 R樹、 R+樹、 CELL樹、四叉樹和八叉樹等索引結構,而在這些結構中KD樹和八叉樹在3D點雲數據組織中應用較為廣泛,PCL對上述兩者進行了實現。
1、八叉樹的描述
八叉樹結構是由 Hunter 博士於1978年首次提出的一種數據模型。八叉樹結構通過對三維空間的幾何實體進行體元剖分,每個體元具有相同的時間和空間復雜度,通過循環遞歸的划分方法對大小為 (2n∗2n∗2n)(2n∗2n∗2n) 的三維空間的幾何對象進行剖分,從而構成一個具有根節點的方向圖。在八叉樹結構中如果被划分的體元具有相同的屬性,則該體元構成一個葉節點;否則繼續對該體元剖分成8個子立方體,依次遞剖分,對於 (2n∗2n∗2n)(2n∗2n∗2n)大小的空間對象,最多剖分 nn 次,如下示意圖所示。
八叉樹結構示意圖
八叉樹(Octree)是一種用於描述三維空間的樹狀數據結構。八叉樹的每個節點表示一個正方體的體積元素,每個節點有八個子節點,這八個子節點所表示的體積元素加在一起就等於父節點的體積。一般中心點作為節點的分叉中心。八叉樹若不為空樹的話,樹中任一節點的子節點恰好只會有八個,或零個,也就是子節點不會有0與8以外的數目。
八叉樹葉子節點代表了分辨率最高的情況。例如分辨率設成0.01m,那么每個葉子就是一個1cm見方的小方塊。
那么,這要用來做什么?想象一個立方體,我們最少可以切成多少個相同等分的小立方體?答案就是8個。再想象我們有一個房間,房間里某個角落藏着一枚金幣,我們想很快的把金幣找出來,聰明的你會怎么做?我們可以把房間當成一個立方體,先切成八個小立方體,然后排除掉沒有放任何東西的小立方體,再把有可能藏金幣的小立方體繼續切八等份….如此下去,平均在Log8(房間內的所有物品數)的時間內就可找到金幣。因此,Octree就是用在3D空間中的場景管理,可以很快地知道物體在3D場景中的位置,或偵測與其它物體是否有碰撞以及是否在可視范圍內。
2、實現Octree的原理
(1). 設定最大遞歸深度
(2). 找出場景的最大尺寸,並以此尺寸建立第一個立方體
(3). 依序將單位元元素丟入能被包含且沒有子節點的立方體
(4). 若沒有達到最大遞歸深度,就進行細分八等份,再將該立方體所裝的單位元元素全部分擔給八個子立方體
(5). 若發現子立方體所分配到的單位元元素數量不為零且跟父立方體是一樣的,則該子立方體停止細分,因為跟據空間分割理論,細分的空間所得到的分配必定較少,若是一樣數目,則再怎么切數目還是一樣,會造成無窮切割的情形。
(6). 重復3,直到達到最大遞歸深度。
3、Octree的存儲結構
本例中Octree的存貯結構用一個有(若干+八)個字段的記錄來表示樹中的每個結點。其中若干字段用來描述該結點的特性(本例中的特性為:節點的值和節點坐標),其余的八個字段用來作為存放指向其八個子結點的指針。此外,還有線性存儲和1托8式存儲。
4、BSP Tree和Octree對比
a) BSP Tree將場景分割為1個面,而Octree分割為3個面。
b) BSP Tree每個節點最多有2個子結點,而Octree最多有8個子結點
5、八叉樹和k-d樹比較
因此BSP Tree可以用在不論幾唯的場景中,而Octree則常用於三維場景
八叉樹算法的算法實現簡單,但大數據量點雲數據下,其使用比較困難的是最小粒度(葉節點)的確定,粒度較大時,有的節點數據量可能仍比較大,后續查詢效率仍比較低,反之,粒度較小,八叉樹的深度增加,需要的內存空間也比較大(每個非葉子節點需要八個指針),效率也降低。而等分的划分依據,使得在數據重心有偏斜的情況下,受划分深度限制,其效率不是太高。
k-d在鄰域查找上比較有優勢,但在大數據量的情況下,若划分粒度較小時,建樹的開銷也較大,但比八叉樹靈活些。在小數據量的情況下,其搜索效率比較高,但在數據量增大的情況下,其效率會有一定的下降,一般是線性上升的規律。
也有將八叉樹和k-d樹結合起來的應用,應用八叉樹進行大粒度的划分和查找,而后使用k-d樹進行細分,效率會有一定的提升,但其搜索效率變化也與數據量的變化有一個線性關系。
6、八叉樹的C++實現
#include <iostream> using namespace std; //定義八叉樹節點類 template<class T> struct OctreeNode { T data; //節點數據 T xmin,xmax; //節點坐標,即六面體個頂點的坐標 T ymin,ymax; T zmin,zmax; OctreeNode <T>*top_left_front,*top_left_back; //該節點的個子結點 OctreeNode <T>*top_right_front,*top_right_back; OctreeNode <T>*bottom_left_front,*bottom_left_back; OctreeNode <T>*bottom_right_front,*bottom_right_back; OctreeNode //節點類 (T nodeValue = T(), T xminValue = T(),T xmaxValue = T(), T yminValue = T(),T ymaxValue = T(), T zminValue = T(),T zmaxValue = T(), OctreeNode<T>*top_left_front_Node = NULL, OctreeNode<T>*top_left_back_Node = NULL, OctreeNode<T>*top_right_front_Node = NULL, OctreeNode<T>*top_right_back_Node = NULL, OctreeNode<T>*bottom_left_front_Node = NULL, OctreeNode<T>*bottom_left_back_Node = NULL, OctreeNode<T>*bottom_right_front_Node = NULL, OctreeNode<T>*bottom_right_back_Node = NULL ) :data(nodeValue), xmin(xminValue),xmax(xmaxValue), ymin(yminValue),ymax(ymaxValue), zmin(zminValue),zmax(zmaxValue), top_left_front(top_left_front_Node), top_left_back(top_left_back_Node), top_right_front(top_right_front_Node), top_right_back(top_right_back_Node), bottom_left_front(bottom_left_front_Node), bottom_left_back(bottom_left_back_Node), bottom_right_front(bottom_right_front_Node), bottom_right_back(bottom_right_back_Node){} }; //創建八叉樹 template <class T> void createOctree(OctreeNode<T> * &root,int maxdepth,double xmin,double xmax,double ymin,double ymax,double zmin,double zmax) { //cout<<"處理中,請稍候……"<<endl; maxdepth=maxdepth-1; //每遞歸一次就將最大遞歸深度-1 if(maxdepth>=0) { root=new OctreeNode<T>(); cout<<"請輸入節點值:"; //root->data =9;//為節點賦值,可以存儲節點信息,如物體可見性。由於是簡單實現八叉樹功能,簡單賦值為9。 cin>>root->data; //為節點賦值 root->xmin=xmin; //為節點坐標賦值 root->xmax=xmax; root->ymin=ymin; root->ymax=ymax; root->zmin=zmin; root->zmax=zmax; double xm=(xmax-xmin)/2;//計算節點個維度上的半邊長 double ym=(ymax-ymin)/2; double zm=(ymax-ymin)/2; //遞歸創建子樹,根據每一個節點所處(是幾號節點)的位置決定其子結點的坐標。 createOctree(root->top_left_front,maxdepth,xmin,xmax-xm,ymax-ym,ymax,zmax-zm,zmax); createOctree(root->top_left_back,maxdepth,xmin,xmax-xm,ymin,ymax-ym,zmax-zm,zmax); createOctree(root->top_right_front,maxdepth,xmax-xm,xmax,ymax-ym,ymax,zmax-zm,zmax); createOctree(root->top_right_back,maxdepth,xmax-xm,xmax,ymin,ymax-ym,zmax-zm,zmax); createOctree(root->bottom_left_front,maxdepth,xmin,xmax-xm,ymax-ym,ymax,zmin,zmax-zm); createOctree(root->bottom_left_back,maxdepth,xmin,xmax-xm,ymin,ymax-ym,zmin,zmax-zm); createOctree(root->bottom_right_front,maxdepth,xmax-xm,xmax,ymax-ym,ymax,zmin,zmax-zm); createOctree(root->bottom_right_back,maxdepth,xmax-xm,xmax,ymin,ymax-ym,zmin,zmax-zm); } } int i=1; //先序遍歷八叉樹 template <class T> void preOrder( OctreeNode<T> * & p) { if(p) { cout<<i<<".當前節點的值為:"<<p->data<<"\n坐標為:"; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; i+=1; cout<<endl; preOrder(p->top_left_front); preOrder(p->top_left_back); preOrder(p->top_right_front); preOrder(p->top_right_back); preOrder(p->bottom_left_front); preOrder(p->bottom_left_back); preOrder(p->bottom_right_front); preOrder(p->bottom_right_back); cout<<endl; } } //求八叉樹的深度 template<class T> int depth(OctreeNode<T> *& p) { if(p == NULL) return -1; int h =depth(p->top_left_front); return h+1; } //計算單位長度,為查找點做准備 int cal(int num) { int result=1; if(1==num) result=1; else { for(int i=1;i<num;i++) result=2*result; } return result; } //查找點 int maxdepth=0; int times=0; static double xmin=0,xmax=0,ymin=0,ymax=0,zmin=0,zmax=0; int tmaxdepth=0; double txm=1,tym=1,tzm=1; template<class T> void find(OctreeNode<T> *& p,double x,double y,double z) { double xm=(p->xmax-p->xmin)/2; double ym=(p->ymax-p->ymin)/2; double zm=(p->ymax-p->ymin)/2; times++; if(x>xmax || x<xmin|| y>ymax || y<ymin || z>zmax || z<zmin) { cout<<"該點不在場景中!"<<endl; return; } if(x<=p->xmin+txm&& x>=p->xmax-txm && y<=p->ymin+tym &&y>=p->ymax-tym && z<=p->zmin+tzm &&z>=p->zmax-tzm ) { cout<<endl<<"找到該點!"<<"該點位於"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<"節點內!"<<endl; cout<<"共經過"<<times<<"次遞歸!"<<endl; } else if(x<(p->xmax-xm) && y<(p->ymax-ym) &&z<(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->bottom_left_back,x,y,z); } else if(x<(p->xmax-xm) && y<(p->ymax-ym) &&z>(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->top_left_back,x,y,z); } else if(x>(p->xmax-xm) && y<(p->ymax-ym) &&z<(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->bottom_right_back,x,y,z); } else if(x>(p->xmax-xm) && y<(p->ymax-ym) &&z>(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->top_right_back,x,y,z); } else if(x<(p->xmax-xm) && y>(p->ymax-ym) &&z<(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->bottom_left_front,x,y,z); } else if(x<(p->xmax-xm) && y>(p->ymax-ym) &&z>(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->top_left_front,x,y,z); } else if(x>(p->xmax-xm) && y>(p->ymax-ym) &&z<(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->bottom_right_front,x,y,z); } else if(x>(p->xmax-xm) && y>(p->ymax-ym) &&z>(p->zmax-zm)) { cout<<"當前經過節點坐標:"<<endl; cout<<"xmin: "<<p->xmin<<" xmax: "<<p->xmax; cout<<"ymin: "<<p->ymin<<" ymax: "<<p->ymax; cout<<"zmin: "<<p->zmin<<" zmax: "<<p->zmax; cout<<endl; find(p->top_right_front,x,y,z); } } //main函數 int main () { OctreeNode<double> *rootNode = NULL; int choiced = 0; while(true) { system("cls"); cout<<"請選擇操作:\n"; cout<<"1.創建八叉樹 2.先序遍歷八叉樹\n"; cout<<"3.查看樹深度 4.查找節點 \n"; cout<<"0.退出\n\n"; cin>>choiced; if(choiced == 0) return 0; else if(choiced == 1) { system("cls"); cout<<"請輸入最大遞歸深度:"<<endl; cin>>maxdepth; cout<<"請輸入外包盒坐標,順序如下:xmin,xmax,ymin,ymax,zmin,zmax"<<endl; cin>>xmin>>xmax>>ymin>>ymax>>zmin>>zmax; if(maxdepth>=0|| xmax>xmin || ymax>ymin || zmax>zmin || xmin>0 || ymin>0||zmin>0) { tmaxdepth=cal(maxdepth); txm=(xmax-xmin)/tmaxdepth; tym=(ymax-ymin)/tmaxdepth; tzm=(zmax-zmin)/tmaxdepth; createOctree(rootNode,maxdepth,xmin,xmax,ymin,ymax,zmin,zmax); } else { cout<<"輸入錯誤!"; return 0; } } else if(choiced == 2) { system("cls"); cout<<"先序遍歷八叉樹結果:/n"; i=1; preOrder(rootNode); cout<<endl; system("pause"); } else if(choiced == 3) { system("cls"); int dep =depth(rootNode); cout<<"此八叉樹的深度為"<<dep+1<<endl; system("pause"); } else if(choiced == 4) { system("cls"); cout<<"請輸入您希望查找的點的坐標,順序如下:x,y,z\n"; double x,y,z; cin>>x>>y>>z; times=0; cout<<endl<<"開始搜尋該點……"<<endl; find(rootNode,x,y,z); system("pause"); } else { system("cls"); cout<<"\n\n錯誤選擇!\n"; system("pause"); } } }
PCL中octree應用
1、八叉樹在PCL中的壓縮例子
// 用於顯示數據 pcl::visualization::CloudViewer viewer; // 用於壓縮、解壓點雲數據 pcl::io::OctreePointCloudCompression<pcl::PointXYZRGBA>* PointCloudEncoder; pcl::io::OctreePointCloudCompression<pcl::PointXYZRGBA>* PointCloudDecoder; // 用於接收到設備數據時的回調響應 void cloud_callback(const pcl::PointCloud<pcl::PointXYZRGBA>::ConstPtr &cloud); // interface->start ();之后,設備獲取一幀數據,就回調一次次函數,不斷刷新壓縮后的點雲 void cloud_callback(const pcl::PointCloud<pcl::PointXYZRGBA>::ConstPtr &cloud) { if (!viewer.wasStopped()) { // 存儲壓縮點雲的字節流對象 std::stringstream compressedData; // 存儲輸出點雲 pcl::PointCloud<pcl::PointXYZRGBA>::Ptr cloudOut(new pcl::PointCloud<pcl::PointXYZRGBA>()); // 壓縮點雲 PointCloudEncoder->encodePointCloud(cloud, compressedData); // 解壓縮點雲 PointCloudDecoder->decodePointCloud(compressedData, cloudOut); // 可視化解壓縮的點雲 viewer.showCloud(cloudOut); } } int main() { // 壓縮選項詳情在: /io/include/pcl/compression/compression_profiles.h pcl::io::compression_Profiles_e compressionProfile = pcl::io::MED_RES_ONLINE_COMPRESSION_WITH_COLOR; // 初始化壓縮和解壓縮對象 其中壓縮對象需要設定壓縮參數選項,解壓縮按照數據源自行判斷 PointCloudEncoder = new pcl::io::OctreePointCloudCompression<pcl::PointXYZRGBA> (compressionProfile, true); PointCloudDecoder = new pcl::io::OctreePointCloudCompression<pcl::PointXYZRGBA> (); //創建從OpenNI獲取點雲的抓取對象 pcl::Grabber* interface = new pcl::OpenNIGrabber (); // 建立回調函數 boost::function<void (const pcl::PointCloud<pcl::PointXYZRGBA>::ConstPtr&)> f = boost::bind (&SimpleOpenNIViewer::cloud_callback, this, _1); //建立回調函數和回調信息的綁定 boost::signals2::connection c = interface->registerCallback (f); // 開始接受點雲的數據流 interface->start (); while (!viewer.wasStopped ()) { sleep (1); } interface->stop (); // 刪除壓縮與解壓縮的實例 delete (PointCloudEncoder); delete (PointCloudDecoder); return 0; }
2、體素內鄰近搜索、K鄰近搜索、半徑內鄰近搜索
void main () { pcl::PintCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<PointXYZ>); cloud->width=1000; cloud->height=1; cloud->pionts.resize(cloud->width*cloud->height); for(int i=0;i<cloud->points.size();i++) { cloud->points[i].x=1024.0f*rand()/(RAND_MAX+1.0f); cloud->points[i].y=1024.0f*rand()/(RAND_MAX+1.0f); cloud->points[i].z=1024.0f*rand()/(RAND_MAX+1.0f); } pcl::octree::OctreePointCloudSearch<pcl::PointXYZ>octree(128.0f)//初始化octree,分辨率為128.0 octree.setInputCloud(cloud); octree.addPointsFromInputCloud();// 根據點雲搭建octree vector<int>piontIndexVec;// 存儲體素鄰近的點的索引 //體素內鄰近搜索 if(octree.voxelSearch(searchPoint,piontIndexVec))//searchPoint為搜索的關鍵點 { for(int i=0;i<piontIndexVec.size();i++) { // 每個點的操作 } } //K鄰近搜索 int k=10; vector<int>knnIndex; vector<float>knnDistance; if(octree.nearestKSearch(searchPoint,k,knnIndex,knnDistance)>0) { // 相關操作 } //半徑內鄰近搜索 float radius=256.0f; vector<int>radiusIndex; vector<float>radiusDistance; if(octree.radiusSearch(searchPoint,radius,radiusIndex,radiusDistance)>0) { // 相關操作 } }
3、無序點雲數據集的空間變化檢測
void main () { pcl::PintCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<PointXYZ>); cloud->width=1000; cloud->height=1; cloud->pionts.resize(cloud->width*cloud->height); for(int i=0;i<cloud->points.size();i++) { cloud->points[i].x=1024.0f*rand()/(RAND_MAX+1.0f); cloud->points[i].y=1024.0f*rand()/(RAND_MAX+1.0f); cloud->points[i].z=1024.0f*rand()/(RAND_MAX+1.0f); } //OctreePointCloudChangeDetector 這個對象允許同時在內存中保管兩個八叉樹,內部使用內存池管理 pcl::octree::OctreePointCloudChangeDetector<pcl::PointXYZ>octree(128.0f)//初始化octree,分辨率為128.0 octree.setInputCloud(cloud); octree.addPointsFromInputCloud();// 根據點雲搭建octree octree.switchBuffers();// 交換緩沖區,使用新的八叉樹結構 pcl::PintCloud<pcl::PointXYZ>::Ptr cloudB(new pcl::PointCloud<PointXYZ>); cloudB->width=1000; cloudB->height=1; cloudB->pionts.resize(cloudB->width*cloudB->height); for(int i=0;i<cloudB->points.size();i++) { cloudB->points[i].x=1024.0f*rand()/(RAND_MAX+1.0f); cloudB->points[i].y=1024.0f*rand()/(RAND_MAX+1.0f); cloudB->points[i].z=1024.0f*rand()/(RAND_MAX+1.0f); } octree.setInputCloud(cloudB); octree.addPointsFromInputCloud();// 根據點雲cloudB搭建octree //getPointIndicesFromNewVoxels()探測兩個點雲中體素間的不同,只能檢測增加的點,不能檢測減少的點 vector<int>newPointIndex; octree.getPointIndicesFromNewVoxels(newPointIndex);//cloudB對比於原始cloud,添加了哪些點,這些點是在cloudB中 for(int i=0;i<newPointIndex.size();i++) { // ~~~相關操作 } }