25 #include <QDomElement> 26 #include <QApplication> 28 #include <QStringList> 42 , mShortCircuit( false )
43 , mInitialisedFlag( false )
44 , mSourceProjection( nullptr )
45 , mDestinationProjection( nullptr )
46 , mSourceDatumTransform( -1 )
47 , mDestinationDatumTransform( -1 )
54 , mShortCircuit( false )
55 , mInitialisedFlag( false )
56 , mSourceProjection( nullptr )
57 , mDestinationProjection( nullptr )
58 , mSourceDatumTransform( -1 )
59 , mDestinationDatumTransform( -1 )
69 , mInitialisedFlag( false )
72 , mSourceProjection( nullptr )
73 , mDestinationProjection( nullptr )
74 , mSourceDatumTransform( -1 )
75 , mDestinationDatumTransform( -1 )
82 , mInitialisedFlag( false )
83 , mSourceProjection( nullptr )
84 , mDestinationProjection( nullptr )
85 , mSourceDatumTransform( -1 )
86 , mDestinationDatumTransform( -1 )
102 , mInitialisedFlag( false )
103 , mSourceProjection( nullptr )
104 , mDestinationProjection( nullptr )
105 , mSourceDatumTransform( -1 )
106 , mDestinationDatumTransform( -1 )
110 mSourceCRS.
createFromId( theSourceSrid, theSourceCRSType );
122 if ( mSourceProjection )
124 pj_free( mSourceProjection );
126 if ( mDestinationProjection )
128 pj_free( mDestinationProjection );
169 mShortCircuit =
true;
181 bool useDefaultDatumTransform = ( mSourceDatumTransform == - 1 && mDestinationDatumTransform == -1 );
185 pj_free( mSourceProjection );
187 if ( !useDefaultDatumTransform )
189 sourceProjString = stripDatumTransform( sourceProjString );
191 if ( mSourceDatumTransform != -1 )
196 pj_free( mDestinationProjection );
198 if ( !useDefaultDatumTransform )
200 destProjString = stripDatumTransform( destProjString );
202 if ( mDestinationDatumTransform != -1 )
207 if ( !useDefaultDatumTransform )
209 addNullGridShifts( sourceProjString, destProjString );
212 mSourceProjection = pj_init_plus( sourceProjString.
toUtf8() );
213 mDestinationProjection = pj_init_plus( destProjString.
toUtf8() );
215 #ifdef COORDINATE_TRANSFORM_VERBOSE 220 mInitialisedFlag =
true;
221 if ( !mDestinationProjection )
223 mInitialisedFlag =
false;
225 if ( !mSourceProjection )
227 mInitialisedFlag =
false;
229 #ifdef COORDINATE_TRANSFORM_VERBOSE 230 if ( mInitialisedFlag )
232 QgsDebugMsg(
"------------------------------------------------------------" );
233 QgsDebugMsg(
"The OGR Coordinate transformation for this layer was set to" );
234 QgsLogger::debug<QgsCoordinateReferenceSystem>(
"Input", mSourceCRS, __FILE__, __FUNCTION__, __LINE__ );
235 QgsLogger::debug<QgsCoordinateReferenceSystem>(
"Output", mDestCRS, __FILE__, __FUNCTION__, __LINE__ );
236 QgsDebugMsg(
"------------------------------------------------------------" );
240 QgsDebugMsg(
"------------------------------------------------------------" );
241 QgsDebugMsg(
"The OGR Coordinate transformation FAILED TO INITIALISE!" );
242 QgsDebugMsg(
"------------------------------------------------------------" );
245 if ( !mInitialisedFlag )
247 QgsDebugMsg(
"Coordinate transformation failed to initialize!" );
254 if ( mSourceCRS == mDestCRS )
258 mShortCircuit =
true;
264 mShortCircuit =
false;
280 if ( mShortCircuit || !mInitialisedFlag )
283 double x = thePoint.
x();
284 double y = thePoint.
y();
317 if ( mShortCircuit || !mInitialisedFlag )
341 #ifdef COORDINATE_TRANSFORM_VERBOSE 354 if ( mShortCircuit || !mInitialisedFlag )
375 double xd =
static_cast< double >( x ), yd = static_cast< double >( y );
384 if ( mShortCircuit || !mInitialisedFlag )
410 if ( mShortCircuit || !mInitialisedFlag )
416 int nVertices = poly.
size();
422 for (
int i = 0; i < nVertices; ++i )
441 for (
int i = 0; i < nVertices; ++i )
453 if ( mShortCircuit || !mInitialisedFlag )
480 if ( mShortCircuit || !mInitialisedFlag )
493 int vectorSize = x.
size();
497 for (
int i = 0; i < vectorSize; ++i )
506 for (
int i = 0; i < vectorSize; ++i )
528 if ( mShortCircuit || !mInitialisedFlag )
541 const int nPoints = 1000;
542 double d = sqrt(( rect.
width() * rect.
height() ) / pow( sqrt( static_cast< double >( nPoints ) ) - 1, 2.0 ) );
543 int nXPoints =
static_cast< int >( ceil( rect.
width() / d ) ) + 1;
544 int nYPoints =
static_cast< int >( ceil( rect.
height() / d ) ) + 1;
560 double dx = rect.
width() /
static_cast< double >( nXPoints - 1 );
561 double dy = rect.
height() /
static_cast< double >( nYPoints - 1 );
565 for (
int i = 0; i < nYPoints ; i++ )
571 for (
int j = 0; j < nXPoints; j++ )
573 x[( i*nXPoints ) + j] = pointX;
574 y[( i*nXPoints ) + j] = pointY;
576 z[( i*nXPoints ) + j] = 0.0;
598 for (
int i = 0; i < nXPoints * nYPoints; i++ )
600 if ( !qIsFinite( x[i] ) || !qIsFinite( y[i] ) )
605 if ( handle180Crossover )
616 if ( handle180Crossover )
637 if ( mShortCircuit || !mInitialisedFlag )
643 "The coordinates can not be reprojected. The CRS is: %1" )
644 .arg( mSourceCRS.
toProj4() ),
tr(
"CRS" ) );
650 "The coordinates can not be reprojected. The CRS is: %1" ).arg( mDestCRS.
toProj4() ),
tr(
"CRS" ) );
654 #ifdef COORDINATE_TRANSFORM_VERBOSE 657 QgsDebugMsg(
QString(
"[[[[[[ Number of points to transform: %1 ]]]]]]" ).arg( numPoints ) );
664 if (( pj_is_latlong( mDestinationProjection ) && ( direction ==
ReverseTransform ) )
665 || ( pj_is_latlong( mSourceProjection ) && ( direction ==
ForwardTransform ) ) )
667 for (
int i = 0; i < numPoints; ++i )
678 projResult = pj_transform( mDestinationProjection, mSourceProjection, numPoints, 0, x, y, z );
682 Q_ASSERT( mSourceProjection );
683 Q_ASSERT( mDestinationProjection );
684 projResult = pj_transform( mSourceProjection, mDestinationProjection, numPoints, 0, x, y, z );
687 if ( projResult != 0 )
692 for (
int i = 0; i < numPoints; ++i )
696 points +=
QString(
"(%1, %2)\n" ).
arg( x[i], 0,
'f' ).
arg( y[i], 0,
'f' );
700 points +=
QString(
"(%1, %2)\n" ).
arg( x[i] * RAD_TO_DEG, 0,
'f' ).
arg( y[i] * RAD_TO_DEG, 0,
'f' );
704 dir = ( direction ==
ForwardTransform ) ?
tr(
"forward transform" ) :
tr(
"inverse transform" );
706 char *srcdef = pj_get_def( mSourceProjection, 0 );
707 char *dstdef = pj_get_def( mDestinationProjection, 0 );
711 "PROJ.4: %3 +to %4\n" 721 QgsDebugMsg(
"Projection failed emitting invalid transform signal: " + msg );
732 if (( pj_is_latlong( mDestinationProjection ) && ( direction ==
ForwardTransform ) )
733 || ( pj_is_latlong( mSourceProjection ) && ( direction ==
ReverseTransform ) ) )
735 for (
int i = 0; i < numPoints; ++i )
742 #ifdef COORDINATE_TRANSFORM_VERBOSE 744 .arg( xorg, 0,
'g', 15 ).arg( yorg, 0,
'g', 15 )
745 .arg( *x, 0,
'g', 15 ).arg( *y, 0,
'g', 15 ) );
752 QgsDebugMsg(
"Reading Coordinate Transform from xml ------------------------!" );
755 mSourceCRS.
readXML( mySrcNode );
758 mDestCRS.
readXML( myDestNode );
776 mSourceCRS.
writeXML( mySourceElement, theDoc );
780 mDestCRS.
writeXML( myDestElement, theDoc );
793 +
"/share/proj/" +
QString( name );
800 void QgsCoordinateTransform::setFinder()
825 return transformations;
831 if ( srcSplit.
size() < 2 || destSplit.
size() < 2 )
833 return transformations;
836 int srcAuthCode = srcSplit.
at( 1 ).toInt();
837 int destAuthCode = destSplit.
at( 1 ).toInt();
839 if ( srcAuthCode == destAuthCode )
841 return transformations;
845 searchDatumTransform(
QString(
"SELECT coord_op_code FROM tbl_datum_transform WHERE source_crs_code=%1 AND target_crs_code=%2 ORDER BY deprecated ASC,preferred DESC" ).arg( srcAuthCode ).arg( destAuthCode ),
848 searchDatumTransform(
QString(
"SELECT coord_op_code FROM tbl_datum_transform WHERE source_crs_code = %1 AND target_crs_code=%2 ORDER BY deprecated ASC,preferred DESC" ).arg( destAuthCode ).arg( srcAuthCode ),
849 reverseDirectTransforms );
851 searchDatumTransform(
QString(
"SELECT coord_op_code FROM tbl_datum_transform WHERE (source_crs_code=%1 AND target_crs_code=%2) OR (source_crs_code=%2 AND target_crs_code=%1) ORDER BY deprecated ASC,preferred DESC" ).arg( srcAuthCode ).arg( 4326 ),
854 searchDatumTransform(
QString(
"SELECT coord_op_code FROM tbl_datum_transform WHERE (source_crs_code=%1 AND target_crs_code=%2) OR (source_crs_code=%2 AND target_crs_code=%1) ORDER BY deprecated ASC,preferred DESC" ).arg( destAuthCode ).arg( 4326 ),
859 for ( ; directIt != directTransforms.
constEnd(); ++directIt )
865 directIt = reverseDirectTransforms.
constBegin();
866 for ( ; directIt != reverseDirectTransforms.
constEnd(); ++directIt )
872 for ( ; srcWgsIt != srcToWgs84.
constEnd(); ++srcWgsIt )
875 for ( ; dstWgsIt != destToWgs84.
constEnd(); ++dstWgsIt )
881 return transformations;
884 QString QgsCoordinateTransform::stripDatumTransform(
const QString& proj4 )
890 for (
int i = 0; i < parameterSplit.
size(); ++i )
892 currentParameter = parameterSplit.
at( i );
893 if ( !currentParameter.
startsWith(
"towgs84", Qt::CaseInsensitive )
894 && !currentParameter.
startsWith(
"nadgrids", Qt::CaseInsensitive ) )
896 newProjString.
append(
'+' );
897 newProjString.
append( currentParameter );
898 newProjString.
append(
' ' );
901 return newProjString;
904 void QgsCoordinateTransform::searchDatumTransform(
const QString& sql,
QList< int >& transforms )
908 if ( openResult != SQLITE_OK )
915 int prepareRes = sqlite3_prepare( db, sql.
toAscii(), sql.
size(), &stmt, nullptr );
916 if ( prepareRes != SQLITE_OK )
918 sqlite3_finalize( stmt );
924 while ( sqlite3_step( stmt ) == SQLITE_ROW )
926 cOpCode =
reinterpret_cast< const char *
>( sqlite3_column_text( stmt, 0 ) );
929 sqlite3_finalize( stmt );
939 if ( openResult != SQLITE_OK )
942 return transformString;
946 QString sql =
QString(
"SELECT coord_op_method_code,p1,p2,p3,p4,p5,p6,p7 FROM tbl_datum_transform WHERE coord_op_code=%1" ).
arg( datumTransform );
947 int prepareRes = sqlite3_prepare( db, sql.
toAscii(), sql.
size(), &stmt, nullptr );
948 if ( prepareRes != SQLITE_OK )
950 sqlite3_finalize( stmt );
952 return transformString;
955 if ( sqlite3_step( stmt ) == SQLITE_ROW )
958 int methodCode = sqlite3_column_int( stmt, 0 );
959 if ( methodCode == 9615 )
961 transformString =
"+nadgrids=" +
QString( reinterpret_cast< const char * >( sqlite3_column_text( stmt, 1 ) ) );
963 else if ( methodCode == 9603 || methodCode == 9606 || methodCode == 9607 )
965 transformString +=
"+towgs84=";
966 double p1 = sqlite3_column_double( stmt, 1 );
967 double p2 = sqlite3_column_double( stmt, 2 );
968 double p3 = sqlite3_column_double( stmt, 3 );
969 double p4 = sqlite3_column_double( stmt, 4 );
970 double p5 = sqlite3_column_double( stmt, 5 );
971 double p6 = sqlite3_column_double( stmt, 6 );
972 double p7 = sqlite3_column_double( stmt, 7 );
973 if ( methodCode == 9603 )
984 sqlite3_finalize( stmt );
986 return transformString;
993 if ( openResult != SQLITE_OK )
1000 QString sql =
QString(
"SELECT epsg_nr,source_crs_code,target_crs_code,remarks,scope,preferred,deprecated FROM tbl_datum_transform WHERE coord_op_code=%1" ).
arg( datumTransform );
1001 int prepareRes = sqlite3_prepare( db, sql.
toAscii(), sql.
size(), &stmt, nullptr );
1002 if ( prepareRes != SQLITE_OK )
1004 sqlite3_finalize( stmt );
1005 sqlite3_close( db );
1009 int srcCrsId, destCrsId;
1010 if ( sqlite3_step( stmt ) != SQLITE_ROW )
1012 sqlite3_finalize( stmt );
1013 sqlite3_close( db );
1017 epsgNr = sqlite3_column_int( stmt, 0 );
1018 srcCrsId = sqlite3_column_int( stmt, 1 );
1019 destCrsId = sqlite3_column_int( stmt, 2 );
1020 remarks =
QString::fromUtf8( reinterpret_cast< const char * >( sqlite3_column_text( stmt, 3 ) ) );
1021 scope =
QString::fromUtf8( reinterpret_cast< const char * >( sqlite3_column_text( stmt, 4 ) ) );
1022 preferred = sqlite3_column_int( stmt, 5 ) != 0;
1023 deprecated = sqlite3_column_int( stmt, 6 ) != 0;
1032 sqlite3_finalize( stmt );
1033 sqlite3_close( db );
1037 void QgsCoordinateTransform::addNullGridShifts(
QString& srcProjString,
QString& destProjString )
1040 if ( mDestinationDatumTransform == -1 && srcProjString.
contains(
"+nadgrids" ) )
1042 destProjString +=
" +nadgrids=@null";
1045 if ( mSourceDatumTransform == -1 && destProjString.
contains(
"+nadgrids" ) )
1047 srcProjString +=
" +nadgrids=@null";
1053 if ( mSourceCRS.
authid().
compare(
"EPSG:3857", Qt::CaseInsensitive ) == 0 && mSourceDatumTransform == -1 )
1055 srcProjString +=
" +nadgrids=@null";
1057 if ( mDestCRS.
authid().
compare(
"EPSG:3857", Qt::CaseInsensitive ) == 0 && mDestinationDatumTransform == -1 )
1059 destProjString +=
" +nadgrids=@null";
A rectangle specified with double values.
QString & append(QChar ch)
const QgsCoordinateReferenceSystem & crsByAuthId(const QString &authid)
Returns the CRS for authid, e.g.
bool isEmpty() const
test if rectangle is empty.
void setMinimal()
Set a rectangle so that min corner is at max and max corner is at min.
bool createFromWkt(const QString &theWkt)
Set up this CRS using a WKT spatial ref sys definition.
QDomNode appendChild(const QDomNode &newChild)
void setXMaximum(double x)
Set the maximum x value.
void push_back(const T &value)
QString attribute(const QString &name, const QString &defValue) const
double yMaximum() const
Get the y maximum value (top side of rectangle)
QString geographicCRSAuthId() const
Returns auth id of related geographic CRS.
QStringList split(const QString &sep, SplitBehavior behavior, Qt::CaseSensitivity cs) const
const T & at(int i) const
bool createFromId(const long theId, CrsType theType=PostgisCrsId)
QString tr(const char *sourceText, const char *disambiguation, int n)
double x() const
Get the x value of the point.
QDomElement toElement() const
bool createFromOgcWmsCrs(QString theCrs)
Set up this CRS from the given OGC CRS.
const char * name() const
QString number(int n, int base)
void combineExtentWith(QgsRectangle *rect)
expand the rectangle so that covers both the original rectangle and the given rectangle ...
bool createFromSrsId(const long theSrsId)
Set up this CRS by fetching the appropriate information from the sqlite backend.
QString fromUtf8(const char *str, int size)
double yMinimum() const
Get the y minimum value (bottom side of rectangle)
double xMaximum() const
Get the x maximum value (right side of rectangle)
#define QgsDebugMsgLevel(str, level)
bool readXML(const QDomNode &theNode)
Restores state from the given Dom node.
void setAttribute(const QString &name, const QString &value)
int toInt(bool *ok, int base) const
bool startsWith(const QString &s, Qt::CaseSensitivity cs) const
static void logMessage(const QString &message, const QString &tag=QString::null, MessageLevel level=WARNING)
add a message to the instance (and create it if necessary)
A class to represent a point.
QDomNode namedItem(const QString &name) const
bool contains(QChar ch, Qt::CaseSensitivity cs) const
bool isValid() const
Returns whether this CRS is correctly initialized and usable.
const T & at(int i) const
bool writeXML(QDomNode &theNode, QDomDocument &theDoc) const
Stores state to the given Dom node in the given document.
Class for storing a coordinate reference system (CRS)
QString authid() const
Returns the authority identifier for the CRS, which includes both the authority (eg EPSG) and the CRS...
double y() const
Get the y value of the point.
static QString srsDbFilePath()
Returns the path to the srs.db file.
QString description() const
Returns the descriptive name of the CRS, eg "WGS 84" or "GDA 94 / Vicgrid94".
Custom exception class for Coordinate Reference System related exceptions.
const_iterator constEnd() const
QDomElement createElement(const QString &tagName)
const_iterator constBegin() const
double width() const
Width of the rectangle.
QString applicationDirPath()
int compare(const QString &other) const
QString toString(bool automaticPrecision=false) const
returns string representation of form xmin,ymin xmax,ymax
QString arg(qlonglong a, int fieldWidth, int base, const QChar &fillChar) const
double xMinimum() const
Get the x minimum value (left side of rectangle)
static QgsCRSCache * instance()
void setXMinimum(double x)
Set the minimum x value.
QByteArray toAscii() const
double height() const
Height of the rectangle.
QString toProj4() const
Returns a Proj4 string representation of this CRS.
QByteArray toUtf8() const