1 Star 1 Fork 0

xiyan-100/shapelib

加入 Gitee
与超过 1200万 开发者一起发现、参与优秀开源项目,私有仓库也完全免费 :)
免费加入
克隆/下载
sbnsearch.c 34.77 KB
一键复制 编辑 原始数据 按行查看 历史
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997
/******************************************************************************
* $Id$
*
* Project: Shapelib
* Purpose: Implementation of search in ESRI SBN spatial index.
* Author: Even Rouault, even dot rouault at mines dash paris dot org
*
******************************************************************************
* Copyright (c) 2012-2014, Even Rouault <even dot rouault at mines-paris dot org>
*
* This software is available under the following "MIT Style" license,
* or at the option of the licensee under the LGPL (see COPYING). This
* option is discussed in more detail in shapelib.html.
*
* --
*
* Permission is hereby granted, free of charge, to any person obtaining a
* copy of this software and associated documentation files (the "Software"),
* to deal in the Software without restriction, including without limitation
* the rights to use, copy, modify, merge, publish, distribute, sublicense,
* and/or sell copies of the Software, and to permit persons to whom the
* Software is furnished to do so, subject to the following conditions:
*
* The above copyright notice and this permission notice shall be included
* in all copies or substantial portions of the Software.
*
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
* OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
* FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
* THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
* DEALINGS IN THE SOFTWARE.
******************************************************************************/
#include "shapefil.h"
#include <math.h>
#include <assert.h>
#include <stdlib.h>
#include <string.h>
SHP_CVSID("$Id$")
#ifndef TRUE
# define TRUE 1
# define FALSE 0
#endif
#ifndef USE_CPL
#if defined(_MSC_VER)
# if _MSC_VER < 1900
# define snprintf _snprintf
# endif
#elif defined(WIN32) || defined(_WIN32)
# ifndef snprintf
# define snprintf _snprintf
# endif
#endif
#endif
#define CACHED_DEPTH_LIMIT 8
#ifdef __cplusplus
#define STATIC_CAST(type,x) static_cast<type>(x)
#define REINTERPRET_CAST(type,x) reinterpret_cast<type>(x)
#define CONST_CAST(type,x) const_cast<type>(x)
#define SHPLIB_NULLPTR nullptr
#else
#define STATIC_CAST(type,x) ((type)(x))
#define REINTERPRET_CAST(type,x) ((type)(x))
#define CONST_CAST(type,x) ((type)(x))
#define SHPLIB_NULLPTR NULL
#endif
#define READ_MSB_INT(ptr) \
STATIC_CAST(int, (((STATIC_CAST(unsigned, (ptr)[0])) << 24) | ((ptr)[1] << 16) | ((ptr)[2] << 8) | (ptr)[3]))
typedef unsigned char uchar;
typedef int coord;
/*typedef uchar coord;*/
typedef struct
{
uchar *pabyShapeDesc; /* Cache of (nShapeCount * 8) bytes of the bins. May be NULL. */
int nBinStart; /* Index of first bin for this node. */
int nShapeCount; /* Number of shapes attached to this node. */
int nBinCount; /* Number of bins for this node. May be 0 if node is empty. */
int nBinOffset; /* Offset in file of the start of the first bin. May be 0 if node is empty. */
int bBBoxInit; /* TRUE if the following bounding box has been computed. */
coord bMinX; /* Bounding box of the shapes directly attached to this node. */
coord bMinY; /* This is *not* the theoretical footprint of the node. */
coord bMaxX;
coord bMaxY;
} SBNNodeDescriptor;
struct SBNSearchInfo
{
SAHooks sHooks;
SAFile fpSBN;
SBNNodeDescriptor *pasNodeDescriptor;
int nShapeCount; /* Total number of shapes */
int nMaxDepth; /* Tree depth */
double dfMinX; /* Bounding box of all shapes */
double dfMaxX;
double dfMinY;
double dfMaxY;
#ifdef DEBUG_IO
int nTotalBytesRead;
#endif
};
typedef struct
{
SBNSearchHandle hSBN;
coord bMinX; /* Search bounding box */
coord bMinY;
coord bMaxX;
coord bMaxY;
int nShapeCount;
int nShapeAlloc;
int *panShapeId; /* 0 based */
uchar abyBinShape[8 * 100];
#ifdef DEBUG_IO
int nBytesRead;
#endif
} SearchStruct;
/************************************************************************/
/* SwapWord() */
/* */
/* Swap a 2, 4 or 8 byte word. */
/************************************************************************/
static void SwapWord( int length, void * wordP )
{
int i;
uchar temp;
for( i=0; i < length/2; i++ )
{
temp = STATIC_CAST(uchar*, wordP)[i];
STATIC_CAST(uchar*, wordP)[i] = STATIC_CAST(uchar*, wordP)[length-i-1];
STATIC_CAST(uchar*, wordP)[length-i-1] = temp;
}
}
/************************************************************************/
/* SBNOpenDiskTree() */
/************************************************************************/
SBNSearchHandle SBNOpenDiskTree( const char* pszSBNFilename,
SAHooks *psHooks )
{
int i;
SBNSearchHandle hSBN;
uchar abyHeader[108];
int nShapeCount;
int nMaxDepth;
int nMaxNodes;
int nNodeDescSize;
int nNodeDescCount;
uchar* pabyData = SHPLIB_NULLPTR;
SBNNodeDescriptor* pasNodeDescriptor = SHPLIB_NULLPTR;
uchar abyBinHeader[8];
int nCurNode;
int nNextNonEmptyNode;
int nExpectedBinId;
int bBigEndian;
/* -------------------------------------------------------------------- */
/* Establish the byte order on this machine. */
/* -------------------------------------------------------------------- */
i = 1;
if( *REINTERPRET_CAST(unsigned char *, &i) == 1 )
bBigEndian = FALSE;
else
bBigEndian = TRUE;
/* -------------------------------------------------------------------- */
/* Initialize the handle structure. */
/* -------------------------------------------------------------------- */
hSBN = STATIC_CAST(SBNSearchHandle,
calloc(sizeof(struct SBNSearchInfo),1));
if (psHooks == SHPLIB_NULLPTR)
SASetupDefaultHooks( &(hSBN->sHooks) );
else
memcpy( &(hSBN->sHooks), psHooks, sizeof(SAHooks) );
hSBN->fpSBN = hSBN->sHooks.FOpen(pszSBNFilename, "rb");
if (hSBN->fpSBN == SHPLIB_NULLPTR)
{
free(hSBN);
return SHPLIB_NULLPTR;
}
/* -------------------------------------------------------------------- */
/* Check file header signature. */
/* -------------------------------------------------------------------- */
if (hSBN->sHooks.FRead(abyHeader, 108, 1, hSBN->fpSBN) != 1 ||
abyHeader[0] != 0 ||
abyHeader[1] != 0 ||
abyHeader[2] != 0x27 ||
(abyHeader[3] != 0x0A && abyHeader[3] != 0x0D) ||
abyHeader[4] != 0xFF ||
abyHeader[5] != 0xFF ||
abyHeader[6] != 0xFE ||
abyHeader[7] != 0x70)
{
hSBN->sHooks.Error( ".sbn file is unreadable, or corrupt." );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* -------------------------------------------------------------------- */
/* Read shapes bounding box. */
/* -------------------------------------------------------------------- */
memcpy(&hSBN->dfMinX, abyHeader + 32, 8);
memcpy(&hSBN->dfMinY, abyHeader + 40, 8);
memcpy(&hSBN->dfMaxX, abyHeader + 48, 8);
memcpy(&hSBN->dfMaxY, abyHeader + 56, 8);
if( !bBigEndian )
{
SwapWord(8, &hSBN->dfMinX);
SwapWord(8, &hSBN->dfMinY);
SwapWord(8, &hSBN->dfMaxX);
SwapWord(8, &hSBN->dfMaxY);
}
if( hSBN->dfMinX > hSBN->dfMaxX ||
hSBN->dfMinY > hSBN->dfMaxY )
{
hSBN->sHooks.Error( "Invalid extent in .sbn file." );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* -------------------------------------------------------------------- */
/* Read and check number of shapes. */
/* -------------------------------------------------------------------- */
nShapeCount = READ_MSB_INT(abyHeader + 28);
hSBN->nShapeCount = nShapeCount;
if (nShapeCount < 0 || nShapeCount > 256000000 )
{
char szErrorMsg[64];
snprintf(szErrorMsg, sizeof(szErrorMsg),
"Invalid shape count in .sbn : %d", nShapeCount );
hSBN->sHooks.Error( szErrorMsg );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* Empty spatial index */
if( nShapeCount == 0 )
{
return hSBN;
}
/* -------------------------------------------------------------------- */
/* Compute tree depth. */
/* It is computed such as in average there are not more than 8 */
/* shapes per node. With a minimum depth of 2, and a maximum of 24 */
/* -------------------------------------------------------------------- */
nMaxDepth = 2;
while( nMaxDepth < 24 && nShapeCount > ((1 << nMaxDepth) - 1) * 8 )
nMaxDepth ++;
hSBN->nMaxDepth = nMaxDepth;
nMaxNodes = (1 << nMaxDepth) - 1;
/* -------------------------------------------------------------------- */
/* Check that the first bin id is 1. */
/* -------------------------------------------------------------------- */
if( READ_MSB_INT(abyHeader + 100) != 1 )
{
hSBN->sHooks.Error( "Unexpected bin id" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* -------------------------------------------------------------------- */
/* Read and check number of node descriptors to be read. */
/* There are at most (2^nMaxDepth) - 1, but all are not necessary */
/* described. Non described nodes are empty. */
/* -------------------------------------------------------------------- */
nNodeDescSize = READ_MSB_INT(abyHeader + 104);
nNodeDescSize *= 2; /* 16-bit words */
/* each bin descriptor is made of 2 ints */
nNodeDescCount = nNodeDescSize / 8;
if ((nNodeDescSize % 8) != 0 ||
nNodeDescCount < 0 || nNodeDescCount > nMaxNodes )
{
char szErrorMsg[64];
snprintf(szErrorMsg, sizeof(szErrorMsg),
"Invalid node descriptor size in .sbn : %d", nNodeDescSize );
hSBN->sHooks.Error( szErrorMsg );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* coverity[tainted_data] */
pabyData = STATIC_CAST(uchar*, malloc( nNodeDescSize ));
pasNodeDescriptor = STATIC_CAST(SBNNodeDescriptor*,
calloc ( nMaxNodes, sizeof(SBNNodeDescriptor) ));
if (pabyData == SHPLIB_NULLPTR || pasNodeDescriptor == SHPLIB_NULLPTR)
{
free(pabyData);
free(pasNodeDescriptor);
hSBN->sHooks.Error( "Out of memory error" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* -------------------------------------------------------------------- */
/* Read node descriptors. */
/* -------------------------------------------------------------------- */
if (hSBN->sHooks.FRead(pabyData, nNodeDescSize, 1,
hSBN->fpSBN) != 1)
{
free(pabyData);
free(pasNodeDescriptor);
hSBN->sHooks.Error( "Cannot read node descriptors" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
hSBN->pasNodeDescriptor = pasNodeDescriptor;
for(i = 0; i < nNodeDescCount; i++)
{
/* -------------------------------------------------------------------- */
/* Each node descriptor contains the index of the first bin that */
/* described it, and the number of shapes in this first bin and */
/* the following bins (in the relevant case). */
/* -------------------------------------------------------------------- */
int nBinStart = READ_MSB_INT(pabyData + 8 * i);
int nNodeShapeCount = READ_MSB_INT(pabyData + 8 * i + 4);
pasNodeDescriptor[i].nBinStart = nBinStart > 0 ? nBinStart : 0;
pasNodeDescriptor[i].nShapeCount = nNodeShapeCount;
if ((nBinStart > 0 && nNodeShapeCount == 0) ||
nNodeShapeCount < 0 || nNodeShapeCount > nShapeCount)
{
hSBN->sHooks.Error( "Inconsistent shape count in bin" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
}
free(pabyData);
/* pabyData = SHPLIB_NULLPTR; */
/* Locate first non-empty node */
nCurNode = 0;
while(nCurNode < nMaxNodes && pasNodeDescriptor[nCurNode].nBinStart <= 0)
nCurNode ++;
if( nCurNode >= nMaxNodes)
{
hSBN->sHooks.Error( "All nodes are empty" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
pasNodeDescriptor[nCurNode].nBinOffset =
STATIC_CAST(int, hSBN->sHooks.FTell(hSBN->fpSBN));
/* Compute the index of the next non empty node. */
nNextNonEmptyNode = nCurNode + 1;
while(nNextNonEmptyNode < nMaxNodes &&
pasNodeDescriptor[nNextNonEmptyNode].nBinStart <= 0)
nNextNonEmptyNode ++;
nExpectedBinId = 1;
/* -------------------------------------------------------------------- */
/* Traverse bins to compute the offset of the first bin of each */
/* node. */
/* Note: we could use the .sbx file to compute the offsets instead.*/
/* -------------------------------------------------------------------- */
while( hSBN->sHooks.FRead(abyBinHeader, 8, 1,
hSBN->fpSBN) == 1 )
{
int nBinId;
int nBinSize;
nExpectedBinId ++;
nBinId = READ_MSB_INT(abyBinHeader);
nBinSize = READ_MSB_INT(abyBinHeader + 4);
nBinSize *= 2; /* 16-bit words */
if( nBinId != nExpectedBinId )
{
hSBN->sHooks.Error( "Unexpected bin id" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
/* Bins are always limited to 100 features */
/* If there are more, then they are located in continuous bins */
if( (nBinSize % 8) != 0 || nBinSize <= 0 || nBinSize > 100 * 8)
{
hSBN->sHooks.Error( "Unexpected bin size" );
SBNCloseDiskTree(hSBN);
return SHPLIB_NULLPTR;
}
if( nNextNonEmptyNode < nMaxNodes &&
nBinId == pasNodeDescriptor[nNextNonEmptyNode].nBinStart )
{
nCurNode = nNextNonEmptyNode;
pasNodeDescriptor[nCurNode].nBinOffset =
STATIC_CAST(int, hSBN->sHooks.FTell(hSBN->fpSBN)) - 8;
/* Compute the index of the next non empty node. */
nNextNonEmptyNode = nCurNode + 1;
while(nNextNonEmptyNode < nMaxNodes &&
pasNodeDescriptor[nNextNonEmptyNode].nBinStart <= 0)
nNextNonEmptyNode ++;
}
pasNodeDescriptor[nCurNode].nBinCount ++;
/* Skip shape description */
hSBN->sHooks.FSeek(hSBN->fpSBN, nBinSize, SEEK_CUR);
}
return hSBN;
}
/***********************************************************************/
/* SBNCloseDiskTree() */
/************************************************************************/
void SBNCloseDiskTree( SBNSearchHandle hSBN )
{
int i;
int nMaxNodes;
if (hSBN == SHPLIB_NULLPTR)
return;
if( hSBN->pasNodeDescriptor != SHPLIB_NULLPTR )
{
nMaxNodes = (1 << hSBN->nMaxDepth) - 1;
for(i = 0; i < nMaxNodes; i++)
{
if( hSBN->pasNodeDescriptor[i].pabyShapeDesc != SHPLIB_NULLPTR )
free(hSBN->pasNodeDescriptor[i].pabyShapeDesc);
}
}
/* printf("hSBN->nTotalBytesRead = %d\n", hSBN->nTotalBytesRead); */
hSBN->sHooks.FClose(hSBN->fpSBN);
free(hSBN->pasNodeDescriptor);
free(hSBN);
}
/************************************************************************/
/* SfRealloc() */
/* */
/* A realloc cover function that will access a NULL pointer as */
/* a valid input. */
/************************************************************************/
static void * SfRealloc( void * pMem, int nNewSize )
{
if( pMem == SHPLIB_NULLPTR )
return malloc(nNewSize);
else
return realloc(pMem,nNewSize);
}
/************************************************************************/
/* SBNAddShapeId() */
/************************************************************************/
static int SBNAddShapeId( SearchStruct* psSearch,
int nShapeId )
{
if (psSearch->nShapeCount == psSearch->nShapeAlloc)
{
int* pNewPtr;
psSearch->nShapeAlloc =
STATIC_CAST(int, ((psSearch->nShapeCount + 100) * 5) / 4);
pNewPtr =
STATIC_CAST(int *, SfRealloc( psSearch->panShapeId,
psSearch->nShapeAlloc * sizeof(int) ));
if( pNewPtr == SHPLIB_NULLPTR )
{
psSearch->hSBN->sHooks.Error( "Out of memory error" );
return FALSE;
}
psSearch->panShapeId = pNewPtr;
}
psSearch->panShapeId[psSearch->nShapeCount] = nShapeId;
psSearch->nShapeCount ++;
return TRUE;
}
/************************************************************************/
/* SBNSearchDiskInternal() */
/************************************************************************/
/* Due to the way integer coordinates are rounded, */
/* we can use a strict intersection test, except when the node */
/* bounding box or the search bounding box is degenerated. */
#define SEARCH_BB_INTERSECTS(_bMinX, _bMinY, _bMaxX, _bMaxY) \
(((bSearchMinX < _bMaxX && bSearchMaxX > _bMinX) || \
((_bMinX == _bMaxX || bSearchMinX == bSearchMaxX) && \
bSearchMinX <= _bMaxX && bSearchMaxX >= _bMinX)) && \
((bSearchMinY < _bMaxY && bSearchMaxY > _bMinY) || \
((_bMinY == _bMaxY || bSearchMinY == bSearchMaxY ) && \
bSearchMinY <= _bMaxY && bSearchMaxY >= _bMinY)))
static int SBNSearchDiskInternal( SearchStruct* psSearch,
int nDepth,
int nNodeId,
coord bNodeMinX,
coord bNodeMinY,
coord bNodeMaxX,
coord bNodeMaxY )
{
SBNSearchHandle hSBN;
SBNNodeDescriptor* psNode;
coord bSearchMinX = psSearch->bMinX;
coord bSearchMinY = psSearch->bMinY;
coord bSearchMaxX = psSearch->bMaxX;
coord bSearchMaxY = psSearch->bMaxY;
hSBN = psSearch->hSBN;
psNode = &(hSBN->pasNodeDescriptor[nNodeId]);
/* -------------------------------------------------------------------- */
/* Check if this node contains shapes that intersect the search */
/* bounding box. */
/* -------------------------------------------------------------------- */
if ( psNode->bBBoxInit &&
!(SEARCH_BB_INTERSECTS(psNode->bMinX, psNode->bMinY,
psNode->bMaxX, psNode->bMaxY)) )
{
/* No intersection, then don't try to read the shapes attached */
/* to this node */
}
/* -------------------------------------------------------------------- */
/* If this node contains shapes that are cached, then read them. */
/* -------------------------------------------------------------------- */
else if (psNode->pabyShapeDesc != SHPLIB_NULLPTR)
{
int j;
uchar* pabyShapeDesc = psNode->pabyShapeDesc;
/* printf("nNodeId = %d, nDepth = %d\n", nNodeId, nDepth); */
for(j = 0; j < psNode->nShapeCount; j++)
{
coord bMinX = pabyShapeDesc[0];
coord bMinY = pabyShapeDesc[1];
coord bMaxX = pabyShapeDesc[2];
coord bMaxY = pabyShapeDesc[3];
if( SEARCH_BB_INTERSECTS(bMinX, bMinY, bMaxX, bMaxY) )
{
int nShapeId;
nShapeId = READ_MSB_INT(pabyShapeDesc + 4);
/* Caution : we count shape id starting from 0, and not 1 */
nShapeId --;
/*printf("shape=%d, minx=%d, miny=%d, maxx=%d, maxy=%d\n",
nShapeId, bMinX, bMinY, bMaxX, bMaxY);*/
if( !SBNAddShapeId( psSearch, nShapeId ) )
return FALSE;
}
pabyShapeDesc += 8;
}
}
/* -------------------------------------------------------------------- */
/* If the node has attached shapes (that are not (yet) cached), */
/* then retrieve them from disk. */
/* -------------------------------------------------------------------- */
else if (psNode->nBinCount > 0)
{
uchar abyBinHeader[8];
int nBinSize, nShapes;
int nShapeCountAcc = 0;
int i, j;
/* printf("nNodeId = %d, nDepth = %d\n", nNodeId, nDepth); */
hSBN->sHooks.FSeek(hSBN->fpSBN, psNode->nBinOffset, SEEK_SET);
if (nDepth < CACHED_DEPTH_LIMIT)
psNode->pabyShapeDesc = STATIC_CAST(uchar*, malloc(psNode->nShapeCount * 8));
for(i = 0; i < psNode->nBinCount; i++)
{
uchar* pabyBinShape;
#ifdef DEBUG_IO
psSearch->nBytesRead += 8;
#endif
if( hSBN->sHooks.FRead(abyBinHeader, 8, 1,
hSBN->fpSBN) != 1)
{
hSBN->sHooks.Error( "I/O error" );
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
return FALSE;
}
if ( READ_MSB_INT(abyBinHeader + 0) != psNode->nBinStart + i )
{
hSBN->sHooks.Error( "Unexpected bin id" );
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
return FALSE;
}
nBinSize = READ_MSB_INT(abyBinHeader + 4);
nBinSize *= 2; /* 16-bit words */
nShapes = nBinSize / 8;
/* Bins are always limited to 100 features */
if( (nBinSize % 8) != 0 || nShapes <= 0 || nShapes > 100)
{
hSBN->sHooks.Error( "Unexpected bin size" );
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
return FALSE;
}
if( nShapeCountAcc + nShapes > psNode->nShapeCount)
{
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
hSBN->sHooks.Error( "Inconsistent shape count for bin" );
return FALSE;
}
if (nDepth < CACHED_DEPTH_LIMIT && psNode->pabyShapeDesc != SHPLIB_NULLPTR)
{
pabyBinShape = psNode->pabyShapeDesc + nShapeCountAcc * 8;
}
else
{
pabyBinShape = psSearch->abyBinShape;
}
#ifdef DEBUG_IO
psSearch->nBytesRead += nBinSize;
#endif
if (hSBN->sHooks.FRead(pabyBinShape, nBinSize, 1,
hSBN->fpSBN) != 1)
{
hSBN->sHooks.Error( "I/O error" );
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
return FALSE;
}
nShapeCountAcc += nShapes;
if (i == 0 && !psNode->bBBoxInit)
{
psNode->bMinX = pabyBinShape[0];
psNode->bMinY = pabyBinShape[1];
psNode->bMaxX = pabyBinShape[2];
psNode->bMaxY = pabyBinShape[3];
}
for(j = 0; j < nShapes; j++)
{
coord bMinX = pabyBinShape[0];
coord bMinY = pabyBinShape[1];
coord bMaxX = pabyBinShape[2];
coord bMaxY = pabyBinShape[3];
if( !psNode->bBBoxInit )
{
#ifdef sanity_checks
/* -------------------------------------------------------------------- */
/* Those tests only check that the shape bounding box in the bin */
/* are consistent (self-consistent and consistent with the node */
/* they are attached to). They are optional however (as far as */
/* the safety of runtime is concerned at least). */
/* -------------------------------------------------------------------- */
if( !(((bMinX < bMaxX ||
(bMinX == 0 && bMaxX == 0) ||
(bMinX == 255 && bMaxX == 255))) &&
((bMinY < bMaxY ||
(bMinY == 0 && bMaxY == 0) ||
(bMinY == 255 && bMaxY == 255)))) ||
bMaxX < bNodeMinX || bMaxY < bNodeMinY ||
bMinX > bNodeMaxX || bMinY > bNodeMaxY )
{
/* printf("shape %d %d %d %d\n", bMinX, bMinY, bMaxX, bMaxY);*/
/* printf("node %d %d %d %d\n", bNodeMinX, bNodeMinY, bNodeMaxX, bNodeMaxY);*/
hSBN->sHooks.Error(
"Invalid shape bounding box in bin" );
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
return FALSE;
}
#endif
if (bMinX < psNode->bMinX) psNode->bMinX = bMinX;
if (bMinY < psNode->bMinY) psNode->bMinY = bMinY;
if (bMaxX > psNode->bMaxX) psNode->bMaxX = bMaxX;
if (bMaxY > psNode->bMaxY) psNode->bMaxY = bMaxY;
}
if( SEARCH_BB_INTERSECTS(bMinX, bMinY, bMaxX, bMaxY) )
{
int nShapeId;
nShapeId = READ_MSB_INT(pabyBinShape + 4);
/* Caution : we count shape id starting from 0, and not 1 */
nShapeId --;
/*printf("shape=%d, minx=%d, miny=%d, maxx=%d, maxy=%d\n",
nShapeId, bMinX, bMinY, bMaxX, bMaxY);*/
if( !SBNAddShapeId( psSearch, nShapeId ) )
return FALSE;
}
pabyBinShape += 8;
}
}
if( nShapeCountAcc != psNode->nShapeCount)
{
free(psNode->pabyShapeDesc);
psNode->pabyShapeDesc = SHPLIB_NULLPTR;
hSBN->sHooks.Error( "Inconsistent shape count for bin" );
return FALSE;
}
psNode->bBBoxInit = TRUE;
}
/* -------------------------------------------------------------------- */
/* Look up in child nodes. */
/* -------------------------------------------------------------------- */
if( nDepth + 1 < hSBN->nMaxDepth )
{
nNodeId = nNodeId * 2 + 1;
if( (nDepth % 2) == 0 ) /* x split */
{
coord bMid = STATIC_CAST(coord, 1 + (STATIC_CAST(int, bNodeMinX) + bNodeMaxX) / 2);
if( bSearchMinX <= bMid - 1 &&
!SBNSearchDiskInternal( psSearch, nDepth + 1, nNodeId + 1,
bNodeMinX, bNodeMinY,
bMid - 1, bNodeMaxY ) )
{
return FALSE;
}
if( bSearchMaxX >= bMid &&
!SBNSearchDiskInternal( psSearch, nDepth + 1, nNodeId,
bMid, bNodeMinY,
bNodeMaxX, bNodeMaxY ) )
{
return FALSE;
}
}
else /* y split */
{
coord bMid = STATIC_CAST(coord, 1 + (STATIC_CAST(int, bNodeMinY) + bNodeMaxY) / 2);
if( bSearchMinY <= bMid - 1 &&
!SBNSearchDiskInternal( psSearch, nDepth + 1, nNodeId + 1,
bNodeMinX, bNodeMinY,
bNodeMaxX, bMid - 1 ) )
{
return FALSE;
}
if( bSearchMaxY >= bMid &&
!SBNSearchDiskInternal( psSearch, nDepth + 1, nNodeId,
bNodeMinX, bMid,
bNodeMaxX, bNodeMaxY ) )
{
return FALSE;
}
}
}
return TRUE;
}
/************************************************************************/
/* compare_ints() */
/************************************************************************/
/* helper for qsort */
static int
compare_ints( const void * a, const void * b)
{
return *REINTERPRET_CAST(const int*, a) - *REINTERPRET_CAST(const int*, b);
}
/************************************************************************/
/* SBNSearchDiskTree() */
/************************************************************************/
int* SBNSearchDiskTree( SBNSearchHandle hSBN,
double *padfBoundsMin, double *padfBoundsMax,
int *pnShapeCount )
{
double dfMinX, dfMinY, dfMaxX, dfMaxY;
double dfDiskXExtent, dfDiskYExtent;
int bMinX, bMinY, bMaxX, bMaxY;
*pnShapeCount = 0;
dfMinX = padfBoundsMin[0];
dfMinY = padfBoundsMin[1];
dfMaxX = padfBoundsMax[0];
dfMaxY = padfBoundsMax[1];
if( dfMinX > dfMaxX || dfMinY > dfMaxY )
return SHPLIB_NULLPTR;
if( dfMaxX < hSBN->dfMinX || dfMaxY < hSBN->dfMinY ||
dfMinX > hSBN->dfMaxX || dfMinY > hSBN->dfMaxY )
return SHPLIB_NULLPTR;
/* -------------------------------------------------------------------- */
/* Compute the search coordinates in [0,255]x[0,255] coord. space */
/* -------------------------------------------------------------------- */
dfDiskXExtent = hSBN->dfMaxX - hSBN->dfMinX;
dfDiskYExtent = hSBN->dfMaxY - hSBN->dfMinY;
if ( dfDiskXExtent == 0.0 )
{
bMinX = 0;
bMaxX = 255;
}
else
{
if( dfMinX < hSBN->dfMinX )
bMinX = 0;
else
{
double dfMinX_255 = (dfMinX - hSBN->dfMinX)
/ dfDiskXExtent * 255.0;
bMinX = STATIC_CAST(int, floor(dfMinX_255 - 0.005));
if( bMinX < 0 ) bMinX = 0;
}
if( dfMaxX > hSBN->dfMaxX )
bMaxX = 255;
else
{
double dfMaxX_255 = (dfMaxX - hSBN->dfMinX)
/ dfDiskXExtent * 255.0;
bMaxX = STATIC_CAST(int, ceil(dfMaxX_255 + 0.005));
if( bMaxX > 255 ) bMaxX = 255;
}
}
if ( dfDiskYExtent == 0.0 )
{
bMinY = 0;
bMaxY = 255;
}
else
{
if( dfMinY < hSBN->dfMinY )
bMinY = 0;
else
{
double dfMinY_255 = (dfMinY - hSBN->dfMinY)
/ dfDiskYExtent * 255.0;
bMinY = STATIC_CAST(int, floor(dfMinY_255 - 0.005));
if( bMinY < 0 ) bMinY = 0;
}
if( dfMaxY > hSBN->dfMaxY )
bMaxY = 255;
else
{
double dfMaxY_255 = (dfMaxY - hSBN->dfMinY)
/ dfDiskYExtent * 255.0;
bMaxY = STATIC_CAST(int, ceil(dfMaxY_255 + 0.005));
if( bMaxY > 255 ) bMaxY = 255;
}
}
/* -------------------------------------------------------------------- */
/* Run the search. */
/* -------------------------------------------------------------------- */
return SBNSearchDiskTreeInteger(hSBN,
bMinX, bMinY, bMaxX, bMaxY,
pnShapeCount);
}
/************************************************************************/
/* SBNSearchDiskTreeInteger() */
/************************************************************************/
int* SBNSearchDiskTreeInteger( SBNSearchHandle hSBN,
int bMinX, int bMinY, int bMaxX, int bMaxY,
int *pnShapeCount )
{
SearchStruct sSearch;
int bRet;
*pnShapeCount = 0;
if( bMinX > bMaxX || bMinY > bMaxY )
return SHPLIB_NULLPTR;
if( bMaxX < 0 || bMaxY < 0 || bMinX > 255 || bMinY > 255 )
return SHPLIB_NULLPTR;
if( hSBN->nShapeCount == 0 )
return SHPLIB_NULLPTR;
/* -------------------------------------------------------------------- */
/* Run the search. */
/* -------------------------------------------------------------------- */
memset( &sSearch, 0, sizeof(sSearch) );
sSearch.hSBN = hSBN;
sSearch.bMinX = STATIC_CAST(coord, bMinX >= 0 ? bMinX : 0);
sSearch.bMinY = STATIC_CAST(coord, bMinY >= 0 ? bMinY : 0);
sSearch.bMaxX = STATIC_CAST(coord, bMaxX <= 255 ? bMaxX : 255);
sSearch.bMaxY = STATIC_CAST(coord, bMaxY <= 255 ? bMaxY : 255);
sSearch.nShapeCount = 0;
sSearch.nShapeAlloc = 0;
sSearch.panShapeId = STATIC_CAST(int*, calloc(1, sizeof(int)));
#ifdef DEBUG_IO
sSearch.nBytesRead = 0;
#endif
bRet = SBNSearchDiskInternal(&sSearch, 0, 0, 0, 0, 255, 255);
#ifdef DEBUG_IO
hSBN->nTotalBytesRead += sSearch.nBytesRead;
/* printf("nBytesRead = %d\n", sSearch.nBytesRead); */
#endif
if( !bRet )
{
free( sSearch.panShapeId );
*pnShapeCount = 0;
return SHPLIB_NULLPTR;
}
*pnShapeCount = sSearch.nShapeCount;
/* -------------------------------------------------------------------- */
/* Sort the id array */
/* -------------------------------------------------------------------- */
qsort(sSearch.panShapeId, *pnShapeCount, sizeof(int), compare_ints);
return sSearch.panShapeId;
}
/************************************************************************/
/* SBNSearchFreeIds() */
/************************************************************************/
void SBNSearchFreeIds( int* panShapeId )
{
free( panShapeId );
}
马建仓 AI 助手
尝试更多
代码解读
代码找茬
代码优化
C
1
https://gitee.com/xiyan100/shapelib.git
git@gitee.com:xiyan100/shapelib.git
xiyan100
shapelib
shapelib
master

搜索帮助

23e8dbc6 1850385 7e0993f3 1850385