Quantum GIS API Documentation
1.7.4
|
00001 /*************************************************************************** 00002 qgsoverlayanalyzer.cpp - QGIS Tools for vector geometry analysis 00003 ------------------- 00004 begin : 8 Nov 2009 00005 copyright : (C) Carson J. Q. Farmer 00006 email : carson.farmer@gmail.com 00007 ***************************************************************************/ 00008 00009 /*************************************************************************** 00010 * * 00011 * This program is free software; you can redistribute it and/or modify * 00012 * it under the terms of the GNU General Public License as published by * 00013 * the Free Software Foundation; either version 2 of the License, or * 00014 * (at your option) any later version. * 00015 * * 00016 ***************************************************************************/ 00017 /* $Id: qgis.h 9774 2008-12-12 05:41:24Z timlinux $ */ 00018 00019 #include "qgsoverlayanalyzer.h" 00020 00021 #include "qgsapplication.h" 00022 #include "qgsfield.h" 00023 #include "qgsfeature.h" 00024 #include "qgslogger.h" 00025 #include "qgscoordinatereferencesystem.h" 00026 #include "qgsvectorfilewriter.h" 00027 #include "qgsvectordataprovider.h" 00028 #include "qgsdistancearea.h" 00029 #include <QProgressDialog> 00030 00031 bool QgsOverlayAnalyzer::intersection( QgsVectorLayer* layerA, QgsVectorLayer* layerB, 00032 const QString& shapefileName, bool onlySelectedFeatures, 00033 QProgressDialog* p ) 00034 { 00035 if ( !layerA && !layerB ) 00036 { 00037 return false; 00038 } 00039 00040 QgsVectorDataProvider* dpA = layerA->dataProvider(); 00041 QgsVectorDataProvider* dpB = layerB->dataProvider(); 00042 if ( !dpA && !dpB ) 00043 { 00044 return false; 00045 } 00046 00047 QGis::WkbType outputType = dpA->geometryType(); 00048 const QgsCoordinateReferenceSystem crs = layerA->crs(); 00049 QgsFieldMap fieldsA = dpA->fields(); 00050 QgsFieldMap fieldsB = dpB->fields(); 00051 combineFieldLists( fieldsA, fieldsB ); 00052 00053 QgsVectorFileWriter vWriter( shapefileName, dpA->encoding(), fieldsA, outputType, &crs ); 00054 QgsFeature currentFeature; 00055 QgsSpatialIndex index; 00056 00057 //take only selection 00058 if ( onlySelectedFeatures ) 00059 { 00060 const QgsFeatureIds selectionB = layerB->selectedFeaturesIds(); 00061 QgsFeatureIds::const_iterator it = selectionB.constBegin(); 00062 for ( ; it != selectionB.constEnd(); ++it ) 00063 { 00064 if ( !layerB->featureAtId( *it, currentFeature, true, true ) ) 00065 { 00066 continue; 00067 } 00068 index.insertFeature( currentFeature ); 00069 } 00070 //use QgsVectorLayer::featureAtId 00071 const QgsFeatureIds selectionA = layerA->selectedFeaturesIds(); 00072 if ( p ) 00073 { 00074 p->setMaximum( selectionA.size() ); 00075 } 00076 QgsFeature currentFeature; 00077 int processedFeatures = 0; 00078 it = selectionA.constBegin(); 00079 for ( ; it != selectionA.constEnd(); ++it ) 00080 { 00081 if ( p ) 00082 { 00083 p->setValue( processedFeatures ); 00084 } 00085 00086 if ( p && p->wasCanceled() ) 00087 { 00088 break; 00089 } 00090 if ( !layerA->featureAtId( *it, currentFeature, true, true ) ) 00091 { 00092 continue; 00093 } 00094 intersectFeature( currentFeature, &vWriter, layerB, &index ); 00095 ++processedFeatures; 00096 } 00097 00098 if ( p ) 00099 { 00100 p->setValue( selectionA.size() ); 00101 } 00102 } 00103 //take all features 00104 else 00105 { 00106 layerB->select( layerB->pendingAllAttributesList(), QgsRectangle(), true, false ); 00107 while ( layerB->nextFeature( currentFeature ) ) 00108 { 00109 index.insertFeature( currentFeature ); 00110 } 00111 QgsFeature currentFeature; 00112 layerA->select( layerA->pendingAllAttributesList(), QgsRectangle(), true, false ); 00113 00114 int featureCount = layerA->featureCount(); 00115 if ( p ) 00116 { 00117 p->setMaximum( featureCount ); 00118 } 00119 int processedFeatures = 0; 00120 00121 while ( layerA->nextFeature( currentFeature ) ) 00122 { 00123 if ( p ) 00124 { 00125 p->setValue( processedFeatures ); 00126 } 00127 if ( p && p->wasCanceled() ) 00128 { 00129 break; 00130 } 00131 intersectFeature( currentFeature, &vWriter, layerB, &index ); 00132 ++processedFeatures; 00133 } 00134 if ( p ) 00135 { 00136 p->setValue( featureCount ); 00137 } 00138 } 00139 return true; 00140 } 00141 00142 void QgsOverlayAnalyzer::intersectFeature( QgsFeature& f, QgsVectorFileWriter* vfw, 00143 QgsVectorLayer* vl, QgsSpatialIndex* index ) 00144 { 00145 QgsGeometry* featureGeometry = f.geometry(); 00146 QgsGeometry* intersectGeometry = 0; 00147 QgsFeature overlayFeature; 00148 00149 if ( !featureGeometry ) 00150 { 00151 return; 00152 } 00153 00154 QList<int> intersects; 00155 intersects = index->intersects( featureGeometry->boundingBox() ); 00156 QList<int>::const_iterator it = intersects.constBegin(); 00157 QgsFeature outFeature; 00158 for ( ; it != intersects.constEnd(); ++it ) 00159 { 00160 if ( !vl->featureAtId( *it, overlayFeature, true, true ) ) 00161 { 00162 continue; 00163 } 00164 00165 if ( featureGeometry->intersects( overlayFeature.geometry() ) ) 00166 { 00167 intersectGeometry = featureGeometry->intersection( overlayFeature.geometry() ); 00168 00169 outFeature.setGeometry( intersectGeometry ); 00170 QgsAttributeMap attributeMapA = f.attributeMap(); 00171 QgsAttributeMap attributeMapB = overlayFeature.attributeMap(); 00172 combineAttributeMaps( attributeMapA, attributeMapB ); 00173 outFeature.setAttributeMap( attributeMapA ); 00174 00175 //add it to vector file writer 00176 if ( vfw ) 00177 { 00178 vfw->addFeature( outFeature ); 00179 } 00180 } 00181 } 00182 } 00183 00184 void QgsOverlayAnalyzer::combineFieldLists( QgsFieldMap& fieldListA, QgsFieldMap fieldListB ) 00185 { 00186 QList<QString> names; 00187 QMap<int, QgsField>::const_iterator j = fieldListA.constBegin(); 00188 while ( j != fieldListA.constEnd() ) 00189 { 00190 names.append( j.value().name() ); 00191 ++j; 00192 } 00193 QMap<int, QgsField>::const_iterator i = fieldListB.constBegin(); 00194 int count = 0; 00195 int fcount = fieldListA.size(); 00196 QgsField field; 00197 while ( i != fieldListB.constEnd() ) 00198 { 00199 field = i.value(); 00200 while ( names.contains( field.name() ) ) 00201 { 00202 QString name = field.name(); 00203 name.append( "_" ).append( QString( count ) ); 00204 field = QgsField( name, field.type() ); 00205 ++count; 00206 } 00207 fieldListA.insert( fcount, field ); 00208 count = 0; 00209 ++fcount; 00210 ++i; 00211 } 00212 } 00213 00214 void QgsOverlayAnalyzer::combineAttributeMaps( QgsAttributeMap& attributeMapA, QgsAttributeMap attributeMapB ) 00215 { 00216 QMap<int, QVariant>::const_iterator i = attributeMapB.constBegin(); 00217 QVariant attribute; 00218 int fcount = attributeMapA.size(); 00219 while ( i != attributeMapB.constEnd() ) 00220 { 00221 attribute = i.value(); 00222 attributeMapA.insert( fcount, attribute ); 00223 ++i; 00224 ++fcount; 00225 } 00226 } 00227