QGIS API Documentation  2.9.0-Master
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
qgspointsample.cpp
Go to the documentation of this file.
1 #include "qgspointsample.h"
2 #include "qgsgeometry.h"
3 #include "qgsspatialindex.h"
4 #include "qgsvectorfilewriter.h"
5 #include "qgsvectorlayer.h"
6 #include <QFile>
7 #include "mersenne-twister.h"
8 
9 
10 QgsPointSample::QgsPointSample( QgsVectorLayer* inputLayer, const QString& outputLayer, QString nPointsAttribute, QString minDistAttribute ): mInputLayer( inputLayer ),
11  mOutputLayer( outputLayer ), mNumberOfPointsAttribute( nPointsAttribute ), mMinDistanceAttribute( minDistAttribute ), mNCreatedPoints( 0 )
12 {
13 }
14 
15 QgsPointSample::QgsPointSample()
16  : mInputLayer( NULL )
17  , mNCreatedPoints( 0 )
18 {
19 }
20 
22 {
23 }
24 
25 int QgsPointSample::createRandomPoints( QProgressDialog* pd )
26 {
27  Q_UNUSED( pd );
28 
29  //create input layer from id (test if polygon, valid)
30  if ( !mInputLayer )
31  {
32  return 1;
33  }
34 
35  if ( mInputLayer->geometryType() != QGis::Polygon )
36  {
37  return 2;
38  }
39 
40  //delete output file if it already exists
41  if ( QFile::exists( mOutputLayer ) )
42  {
44  }
45 
46  //create vector file writer
47  QgsFields outputFields;
48  outputFields.append( QgsField( "id", QVariant::Int ) );
49  outputFields.append( QgsField( "station_id", QVariant::Int ) );
50  outputFields.append( QgsField( "stratum_id", QVariant::Int ) );
51  QgsVectorFileWriter writer( mOutputLayer, "UTF-8",
52  outputFields,
54  &( mInputLayer->crs() ) );
55 
56  //check if creation of output layer successfull
57  if ( writer.hasError() != QgsVectorFileWriter::NoError )
58  {
59  return 3;
60  }
61 
62  //init random number generator
63  mt_srand( QTime::currentTime().msec() );
64  QgsFeature fet;
65  int nPoints = 0;
66  double minDistance = 0;
67  mNCreatedPoints = 0;
68 
69  QgsFeatureIterator fIt = mInputLayer->getFeatures( QgsFeatureRequest().setSubsetOfAttributes(
70  QStringList() << mNumberOfPointsAttribute << mMinDistanceAttribute, mInputLayer->pendingFields() ) );
71  while ( fIt.nextFeature( fet ) )
72  {
73  nPoints = fet.attribute( mNumberOfPointsAttribute ).toInt();
74  if ( !mMinDistanceAttribute.isEmpty() )
75  {
76  minDistance = fet.attribute( mMinDistanceAttribute ).toDouble();
77  }
78  addSamplePoints( fet, writer, nPoints, minDistance );
79  }
80 
81  return 0;
82 }
83 
84 void QgsPointSample::addSamplePoints( QgsFeature& inputFeature, QgsVectorFileWriter& writer, int nPoints, double minDistance )
85 {
86  QgsGeometry* geom = inputFeature.geometry();
87  if ( !geom )
88  {
89  return;
90  }
91 
92  QgsRectangle geomRect = geom->boundingBox();
93  if ( geomRect.isEmpty() )
94  {
95  return;
96  }
97 
98  QgsSpatialIndex sIndex; //to check minimum distance
99  QMap< QgsFeatureId, QgsPoint > pointMapForFeature;
100 
101  int nIterations = 0;
102  int maxIterations = nPoints * 200;
103  int points = 0;
104 
105  double randX = 0;
106  double randY = 0;
107 
108  while ( nIterations < maxIterations && points < nPoints )
109  {
110  randX = (( double )mt_rand() / MD_RAND_MAX ) * geomRect.width() + geomRect.xMinimum();
111  randY = (( double )mt_rand() / MD_RAND_MAX ) * geomRect.height() + geomRect.yMinimum();
112  QgsPoint randPoint( randX, randY );
113  QgsGeometry* ptGeom = QgsGeometry::fromPoint( randPoint );
114  if ( ptGeom->within( geom ) && checkMinDistance( randPoint, sIndex, minDistance, pointMapForFeature ) )
115  {
116  //add feature to writer
117  QgsFeature f( mNCreatedPoints );
118  f.setAttribute( "id", mNCreatedPoints + 1 );
119  f.setAttribute( "station_id", points + 1 );
120  f.setAttribute( "stratum_id", inputFeature.id() );
121  f.setGeometry( ptGeom );
122  writer.addFeature( f );
123  sIndex.insertFeature( f );
124  pointMapForFeature.insert( mNCreatedPoints, randPoint );
125  ++points;
126  ++mNCreatedPoints;
127  }
128  else
129  {
130  delete ptGeom;
131  }
132  ++nIterations;
133  }
134 }
135 
136 bool QgsPointSample::checkMinDistance( QgsPoint& pt, QgsSpatialIndex& index, double minDistance, QMap< QgsFeatureId, QgsPoint >& pointMap )
137 {
138  if ( minDistance <= 0 )
139  {
140  return true;
141  }
142 
143  QList<QgsFeatureId> neighborList = index.nearestNeighbor( pt, 1 );
144  if ( neighborList.isEmpty() )
145  {
146  return true;
147  }
148 
149  QMap< QgsFeatureId, QgsPoint >::const_iterator it = pointMap.find( neighborList[0] );
150  if ( it == pointMap.constEnd() ) //should not happen
151  {
152  return true;
153  }
154 
155  QgsPoint neighborPt = it.value();
156  if ( neighborPt.sqrDist( pt ) < ( minDistance * minDistance ) )
157  {
158  return false;
159  }
160  return true;
161 }
162 
163 
164 
165 
QgsFeatureId id() const
Get the feature id for this feature.
Definition: qgsfeature.cpp:100
Wrapper for iterator of features from vector data provider or vector layer.
static unsigned index
A rectangle specified with double values.
Definition: qgsrectangle.h:35
bool isEmpty() const
test if rectangle is empty.
int createRandomPoints(QProgressDialog *pd)
Starts calculation of random points.
QgsFeatureIterator getFeatures(const QgsFeatureRequest &request=QgsFeatureRequest())
Query the provider for features specified in request.
QgsGeometry * geometry() const
Get the geometry object associated with this feature.
Definition: qgsfeature.cpp:112
static bool deleteShapeFile(QString theFileName)
Delete a shapefile (and its accompanying shx / dbf / prf)
Container of fields for a vector layer.
Definition: qgsfield.h:172
A convenience class for writing vector files to disk.
The feature class encapsulates a single feature including its id, geometry and a list of field/values...
Definition: qgsfeature.h:113
double sqrDist(double x, double y) const
Returns the squared distance between this point and x,y.
Definition: qgspoint.cpp:333
#define MD_RAND_MAX
int mt_rand()
QgsPointSample(QgsVectorLayer *inputLayer, const QString &outputLayer, QString nPointsAttribute, QString minDistAttribute=QString())
double yMinimum() const
Get the y minimum value (bottom side of rectangle)
Definition: qgsrectangle.h:193
This class wraps a request for features to a vector layer (or directly its vector data provider)...
QList< QgsFeatureId > nearestNeighbor(QgsPoint point, int neighbors) const
returns nearest neighbors (their count is specified by second parameter)
void mt_srand(unsigned value)
bool addFeature(QgsFeature &feature, QgsFeatureRendererV2 *renderer=0, QGis::UnitType outputUnit=QGis::Meters)
add feature to the currently opened shapefile
bool append(const QgsField &field, FieldOrigin origin=OriginProvider, int originIndex=-1)
Append a field. The field must have unique name, otherwise it is rejected (returns false) ...
Definition: qgsfield.cpp:177
QGis::GeometryType geometryType() const
Returns point, line or polygon.
Encapsulate a field in an attribute table or data source.
Definition: qgsfield.h:33
A class to represent a point.
Definition: qgspoint.h:63
static QgsGeometry * fromPoint(const QgsPoint &point)
construct geometry from a point
QgsRectangle boundingBox()
Returns the bounding box of this feature.
bool insertFeature(const QgsFeature &f)
add feature to index
WriterError hasError()
checks whether there were any errors in constructor
QVariant attribute(const QString &name) const
Lookup attribute value from attribute name.
Definition: qgsfeature.cpp:230
bool within(const QgsGeometry *geometry) const
Test for if geometry is within another (uses GEOS)
const QgsCoordinateReferenceSystem & crs() const
Returns layer's spatial reference system.
const QgsFields & pendingFields() const
returns field list in the to-be-committed state
bool nextFeature(QgsFeature &f)
double width() const
Width of the rectangle.
Definition: qgsrectangle.h:198
Represents a vector layer which manages a vector based data sets.
double xMinimum() const
Get the x minimum value (left side of rectangle)
Definition: qgsrectangle.h:183
double height() const
Height of the rectangle.
Definition: qgsrectangle.h:203