做GIS开发的朋友可能对shp并不陌生,但是看到CSDN网友不断提问关于shp文件的一些问题,利用闲暇我对shp文件的一些知识加以总结,共享CSDN网友。
首先了解一下shp文件的一些简单知识
Shapefile文件是美国环境系统研究所(ESRI)所研制的GIS文件系统格式文件,是工业标准的矢量数据文件。 Shapefile将空间特征表中的非拓扑几何对象和属性信息存储在数据集中,特征表中的几何对象存为以坐标点集表示的图形文件—SHP文件,Shapefile文件并不含拓扑(Topological)数据结构。一个Shape文件包括三个文件:一个主文件(*.shp),一个索引文件(*.shx),和一个dBASE(*.dbf)表。主文件是一个直接存取,变长度记录的文件,其中每个记录描述构成一个地理特征(Feature)的所有vertices坐标值。在索引文件中,每条记录包含对应主文件记录距离主文件头开始的偏移量,dBASE表包含SHP文件中每一个Feature的特征属性,表中几何记录和属性数据之间的一一对应关系是基于记录数目的ID。在dBASE文件中的属性记录必须和主文件中的记录顺序是相同的。图形数据和属性数据通过索引号建立一一对应的关系。
Shapefile中坐标文件(.shp)由固定长度的文件头和接着的变长度空间数据记录组成。文件头由100字节的说明信息组成的,主要说明文件的长度、Shape类型、整个Shape图层的范围等等,这些信息构成了空间数据的元数据。在导入空间数据时首先要读入文件头获取Shape文件的基本信息,并以此信息为基础建立相应的元数据表。而变长度空间数据记录是由固定长度的记录头和变长度记录内容组成,其记录结构基本类似,每条记录都有记录头和记录内容组成(空间坐标对)。记录头的内容包括记录号(Record Number)和坐标记录长度(Content Length)两个记录项,Shapefile文件中的记录号都是从1开始的,坐标记录长度是按16位字来衡量的。记录内容包括目标的几何类型(ShapeType)和具体的坐标记录(X,Y),记录内容因要素几何类型的不同,其具体的内容和格式都有所不同。对于具体的记录主要包括空Shape记录,点记录,线记录和多边形记录。
属性文件(.dbf)用于记录属性信息。它是一个标准的DBF文件,也是由头文件和实体信息两部分构成。其中文件头部分的长度是不定长的,它主要对DBF文件作了一些总体说明,其中最主要的是对这个DBF文件的记录项的信息进行了详细的描述,比如对每个记录项的名称,数据类型,长度等信息都有具体的说明。属性文件的实体信息部分就是一条条属性记录,每条记录都是由若干个记录项构成,因此只要依次循环读取每条记录就可以了。
索引文件(.shx)主要包含坐标文件的索引信息,文件中每个记录包含对应的坐标文件记录距离坐标文件的文件头的偏移量。通过索引文件可以很方便地在坐标文件中定位到指定目标地坐标信息。索引文件也是由文件头和实体信息两部分构成的,其中文件头部分是一个长度固定(100 bytes)的记录段,其内容与坐标文件的文件头基本一致。它的实体信息以记录为基本单位,每一条记录包括偏移量(Offset)和记录段长度(Content Length)两个记录项。
接下来我们再看一下代码(c++): SHP文件的读取
- FeatureClass* CGISMapDoc::ImportShapeFileData( FILE* fpShp, FILE* fpDbf )
- {
-
- int fileCode = -1;
- int fileLength = -1;
- int version = -1;
- int shapeType = -1;
- fread(&fileCode , sizeof(int) , 1 , fpShp) ;
- fileCode = ReverseBytes(fileCode) ;
-
- if (fileCode != 9994)
- {
- CString strTemp ;
- strTemp.Format(" WARNING filecode %d ", fileCode );
- AfxMessageBox(strTemp);
- }
-
- for( int i = 0 ; i < 5 ; i ++ )
- fread(&fileCode , sizeof(int) , 1 , fpShp) ;
-
- fread(&fileLength , sizeof(int) , 1 , fpShp) ;
- fileLength = ReverseBytes(fileLength) ;
-
- fread(&version , sizeof(int) , 1 , fpShp) ;
- fread(&shapeType , sizeof(int) , 1 , fpShp) ;
-
- double tempOriginX , tempOriginY ;
- fread( &tempOriginX , sizeof(double) , 1 , fpShp ) ;
- fread( &tempOriginY , sizeof(double) , 1 , fpShp ) ;
-
- double xMaxLayer , yMaxLayer ;
- fread( &xMaxLayer , sizeof(double) , 1 , fpShp ) ;
- fread( &yMaxLayer , sizeof(double) , 1 , fpShp ) ;
-
- double* skip = new double[4] ;
- fread( skip , sizeof(double) , 4 , fpShp ) ;
- delete []skip ;
- skip = 0 ;
-
-
- int uniqueID = this->m_pDataSource->GetUniqueID() ;
- FeatureClass* pShpDataSet = 0 ;
-
- switch( shapeType )
- {
- case 1 :
- pShpDataSet = (FeatureClass*)&(m_pDataSource->CreateDataSet(uniqueID , POINTDATASET , layerName)) ;
- break ;
- case 3 :
- case 23 :
- pShpDataSet = (FeatureClass*)&(m_pDataSource->CreateDataSet(uniqueID , LINEDATASET , layerName)) ;
- break ;
- case 5 :
- pShpDataSet = (FeatureClass*)&(m_pDataSource->CreateDataSet(uniqueID , POLYGONDATASET , layerName)) ;
- break ;
- }
-
- if ( pShpDataSet == 0 ) return 0;
-
-
- struct DBFHeader
- {
- char m_nValid;
- char m_aDate[3];
- char m_nNumRecords[4];
- char m_nHeaderBytes[2];
- char m_nRecordBytes[2];
- char m_nReserved1[3];
- char m_nReserved2[13];
- char m_nReserved3[4];
- }dbfheader;
-
- struct DBFFIELDDescriptor
- {
- char m_sName[10];
- char m_nType;
- char m_nAddress[4];
- char m_nFieldLength;
- char m_nFieldDecimal;
- char m_nReserved1[2];
- char m_nWorkArea;
- char m_nReserved2[2];
- char m_nSetFieldsFlag;
- char m_nReserved3[8];
- };
-
- fread(&dbfheader,sizeof(DBFHeader),1,fpDbf);
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- Fields& fields = pShpDataSet->GetFields();
- DBFFIELDDescriptor field ;
- BYTE endByte = ' ';
- char fieldName[12];
- int fieldDecimal, fieldLen, everyRecordLen = 0 ;
- while ( !feof(fpDbf) )
- {
- fread(&endByte,sizeof(BYTE),1,fpDbf);
- if ( endByte == 0x0D) break ;
- fread(&field,sizeof(DBFFIELDDescriptor),1,fpDbf);
-
- fieldName[0] = endByte;
- for ( i = 0; i < 10; i ++ )
- fieldName[i+1] = field.m_sName[i];
- fieldName[11] = '/0';
-
- fieldDecimal = field.m_nFieldDecimal;
- fieldLen = field.m_nFieldLength;
- switch( field.m_nType )
- {
- case 'C':
- fields.AddField(fieldName,fieldName,FIELD_STRING,fieldLen);
- break;
- case 'F':
- fields.AddField(fieldName,fieldName,FIELD_DOUBLE,fieldLen);
- break;
- case 'N':
- {
- if ( fieldDecimal == 0 )
- fields.AddField(fieldName,fieldName,FIELD_INT,fieldLen);
- else fields.AddField(fieldName,fieldName,FIELD_DOUBLE,fieldLen);
- }
- break;
- }
- everyRecordLen += fieldLen ;
- }
-
-
- while( !feof(fpShp) )
- {
-
- int recordNumber = -1 ;
- int contentLength = -1 ;
- fread(&recordNumber , sizeof(int) , 1 , fpShp) ;
- fread(&contentLength , sizeof(int) , 1 , fpShp) ;
- recordNumber = ReverseBytes(recordNumber) ;
- contentLength = ReverseBytes(contentLength) ;
-
-
- switch( shapeType )
- {
- case 1:
-
- {
- Fields &featureFields = pShpDataSet->GetFields();
- Feature *pFeature = new Feature(recordNumber , 1 , &featureFields) ;
-
- int pointShapeType ;
- fread(&pointShapeType , sizeof(int) , 1 , fpShp) ;
- double xValue , yValue ;
- fread(&xValue , sizeof(double) , 1 , fpShp) ;
- fread(&yValue , sizeof(double) , 1 , fpShp) ;
-
- GeoPoint *pNewGeoPoint = new GeoPoint( xValue , yValue ) ;
- pFeature->SetBound(xValue , yValue , 0 , 0 ) ;
- pFeature->SetGeometry(pNewGeoPoint) ;
- this->LoadAttributeData(pFeature,fpDbf,everyRecordLen);
- pShpDataSet->AddRow(pFeature) ;
- }
- break ;
-
-
- case 3:
-
- {
- Fields &featureFields = pShpDataSet->GetFields();
- Feature *pFeature = new Feature(recordNumber , 1 , &featureFields) ;
-
- int arcShapeType ;
- fread(&arcShapeType , sizeof(int) , 1 , fpShp) ;
-
- double objMinX , objMinY , objMaxX , objMaxY ;
- fread(&objMinX , sizeof(double) , 1 , fpShp) ;
- fread(&objMinY , sizeof(double) , 1 , fpShp) ;
- fread(&objMaxX , sizeof(double) , 1 , fpShp) ;
- fread(&objMaxY , sizeof(double) , 1 , fpShp) ;
-
- GeoPolyline *pNewGeoLine = new GeoPolyline();
- double width = objMaxX - objMinX ;
- double height = objMaxY - objMinY ;
- pFeature->SetBound(objMinX , objMinY , width , height) ;
-
- int numParts , numPoints ;
- fread(&numParts , sizeof(int) , 1 , fpShp) ;
- fread(&numPoints , sizeof(int) , 1 , fpShp) ;
-
- int* startOfPart = new int[numParts] ;
- for( int i = 0 ; i < numParts ; i++ )
- {
- int indexFirstPoint ;
- fread(&indexFirstPoint , sizeof(int) , 1 , fpShp) ;
- startOfPart[i] = indexFirstPoint ;
- }
-
-
- pNewGeoLine->SetPointsCount( numParts ) ;
-
- for( i = 0 ; i < numParts ; i++ )
- {
- GeoPoints& points = pNewGeoLine->GetPoints(i) ;
- int curPosIndex = startOfPart[i] ;
- int nextPosIndex = 0 ;
- int curPointCount = 0 ;
- if( i == numParts - 1 )
- curPointCount = numPoints - curPosIndex ;
- else
- {
- nextPosIndex = startOfPart[i + 1] ;
- curPointCount = nextPosIndex - curPosIndex ;
- }
- points.SetPointCount( curPointCount ) ;
-
- for( int iteratorPoint = 0 ; iteratorPoint < curPointCount ; iteratorPoint ++ )
- {
- double x , y ;
- fread(&x , sizeof(double) , 1 , fpShp) ;
- fread(&y , sizeof(double) , 1 , fpShp) ;
- GeoPoint newVertex(x, y);
- points.SetPoint(iteratorPoint, newVertex);
- }
- }
- delete []startOfPart ;
- startOfPart = 0 ;
- pFeature->SetGeometry(pNewGeoLine) ;
- this->LoadAttributeData(pFeature,fpDbf,everyRecordLen);
- pShpDataSet->AddRow(pFeature) ;
- }
- break ;
-
-
- case 5:
-
- {
- Fields &featureFields = pShpDataSet->GetFields();
- Feature *pFeature = new Feature(recordNumber , 1 , &featureFields) ;
- int polygonShapeType ;
- fread(&polygonShapeType , sizeof(int) , 1 ,fpShp) ;
- if( polygonShapeType != 5 )
- AfxMessageBox( "Error: Attempt to load non polygon shape as polygon." ) ;
-
- double objMinX , objMinY , objMaxX , objMaxY ;
- fread(&objMinX , sizeof(double) , 1 , fpShp) ;
- fread(&objMinY , sizeof(double) , 1 , fpShp) ;
- fread(&objMaxX , sizeof(double) , 1 , fpShp) ;
- fread(&objMaxY , sizeof(double) , 1 , fpShp) ;
-
- GeoPolygon *pNewGeoPolygon = new GeoPolygon();
- double width = objMaxX - objMinX ;
- double height = objMaxY - objMinY ;
- pFeature->SetBound(objMinX , objMinY , width , height) ;
-
- int numParts , numPoints ;
- fread(&numParts , sizeof(int) , 1 , fpShp) ;
- fread(&numPoints , sizeof(int) , 1 , fpShp) ;
-
- int* startOfPart = new int[numParts] ;
- for( int i = 0 ; i < numParts ; i++ )
- {
- int indexFirstPoint ;
- fread(&indexFirstPoint , sizeof(int) , 1 , fpShp) ;
- startOfPart[i] = indexFirstPoint ;
- }
-
-
- pNewGeoPolygon->SetPointsCount( numParts ) ;
-
- for( i = 0 ; i < numParts ; i++ )
- {
- GeoPoints& points = pNewGeoPolygon->GetPoints(i) ;
- int curPosIndex = startOfPart[i] ;
- int nextPosIndex = 0 ;
- int curPointCount = 0 ;
- if( i == numParts - 1 )
- curPointCount = numPoints - curPosIndex ;
- else
- {
- nextPosIndex = startOfPart[i + 1];
- curPointCount = nextPosIndex - curPosIndex ;
- }
- points.SetPointCount( curPointCount ) ;
-
- for( int iteratorPoint = 0 ; iteratorPoint < curPointCount ; iteratorPoint ++ )
- {
- double x , y ;
- fread(&x , sizeof(double) , 1 , fpShp) ;
- fread(&y , sizeof(double) , 1 , fpShp) ;
- GeoPoint newVertex(x, y);
- points.SetPoint(iteratorPoint, newVertex);
- }
- }
- delete []startOfPart ;
- startOfPart = 0 ;
- pFeature->SetGeometry(pNewGeoPolygon) ;
- this->LoadAttributeData(pFeature,fpDbf,everyRecordLen);
- pShpDataSet->AddRow(pFeature) ;
- }
- break ;
-
-
- case 23:
-
- {
- Fields &featureFields = pShpDataSet->GetFields();
- Feature *pFeature = new Feature(recordNumber , 1 , &featureFields) ;
- int arcMShapeType ;
- fread(&arcMShapeType , sizeof(int) , 1 , fpShp) ;
-
- double objMinX , objMinY , objMaxX , objMaxY ;
- fread(&objMinX , sizeof(double) , 1 , fpShp) ;
- fread(&objMinY , sizeof(double) , 1 , fpShp) ;
- fread(&objMaxX , sizeof(double) , 1 , fpShp) ;
- fread(&objMaxY , sizeof(double) , 1 , fpShp) ;
-
- GeoPolyline *pNewGeoLine = new GeoPolyline();
- double width = objMaxX - objMinX ;
- double height = objMaxY - objMinY ;
- pFeature->SetBound(objMinX , objMinY , width , height) ;
-
- int numParts , numPoints ;
- fread(&numParts , sizeof(int) , 1 , fpShp) ;
- fread(&numPoints , sizeof(int) , 1 , fpShp) ;
-
- int* startOfPart = new int[numParts] ;
- for( int i = 0 ; i < numParts ; i++ )
- {
- int indexFirstPoint ;
- fread(&indexFirstPoint , sizeof(int) , 1 , fpShp) ;
- startOfPart[i] = indexFirstPoint ;
- }
-
-
- pNewGeoLine->SetPointsCount( numParts ) ;
-
- for( i = 0 ; i < numParts ; i++ )
- {
- GeoPoints& points = pNewGeoLine->GetPoints(i) ;
- int curPosIndex = startOfPart[i] ;
- int nextPosIndex = 0 ;
- int curPointCount = 0 ;
- if( i == numParts - 1 )
- curPointCount = numPoints - curPosIndex ;
- else
- {
- nextPosIndex = startOfPart[i + 1] ;
- curPointCount = nextPosIndex - curPosIndex ;
- }
- points.SetPointCount( curPointCount ) ;
-
- for( int iteratorPoint = 0 ; iteratorPoint < curPointCount ; iteratorPoint ++ )
- {
- double x , y ;
- fread(&x , sizeof(double) , 1 , fpShp) ;
- fread(&y , sizeof(double) , 1 , fpShp) ;
- GeoPoint newVertex(x, y);
- points.SetPoint(iteratorPoint, newVertex);
- }
- }
- delete []startOfPart ;
- startOfPart = 0 ;
-
- double* value = new double[2 + numPoints] ;
- fread( value , sizeof(double) , 2+numPoints, fpShp) ;
- delete []value ;
- value = 0 ;
-
- pFeature->SetGeometry(pNewGeoLine) ;
- this->LoadAttributeData(pFeature,fpDbf,everyRecordLen);
- pShpDataSet->AddRow(pFeature);
- }
- break ;
-
- }
- }
- return pShpDataSet ;
- }
|