~ubuntu-branches/debian/sid/gdal/sid

« back to all changes in this revision

Viewing changes to ogr/ogrct.cpp

  • Committer: Package Import Robot
  • Author(s): Francesco Paolo Lovergine
  • Date: 2012-05-07 15:04:42 UTC
  • mfrom: (5.5.16 experimental)
  • Revision ID: package-import@ubuntu.com-20120507150442-2eks97loeh6rq005
Tags: 1.9.0-1
* Ready for sid, starting transition.
* All symfiles updated to latest builds.
* Added dh_numpy call in debian/rules to depend on numpy ABI.
* Policy bumped to 3.9.3, no changes required.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
1
/******************************************************************************
2
 
 * $Id: ogrct.cpp 18520 2010-01-11 03:59:14Z warmerdam $
 
2
 * $Id: ogrct.cpp 20079 2010-07-16 21:28:02Z rouault $
3
3
 *
4
4
 * Project:  OpenGIS Simple Features Reference Implementation
5
5
 * Purpose:  The OGRSCoordinateTransformation class.
38
38
#include "proj_api.h"
39
39
#endif
40
40
 
41
 
CPL_CVSID("$Id: ogrct.cpp 18520 2010-01-11 03:59:14Z warmerdam $");
 
41
CPL_CVSID("$Id: ogrct.cpp 20079 2010-07-16 21:28:02Z rouault $");
42
42
 
43
43
/* ==================================================================== */
44
44
/*      PROJ.4 interface stuff.                                         */
47
47
typedef struct { double u, v; } projUV;
48
48
 
49
49
#define projPJ void *
50
 
 
 
50
#define projCtx void *
51
51
#define RAD_TO_DEG      57.29577951308232
52
52
#define DEG_TO_RAD      .0174532925199432958
53
53
 
 
54
#else
 
55
 
 
56
#if PJ_VERSION < 480
 
57
#define projCtx void *
 
58
#endif
 
59
 
54
60
#endif
55
61
 
56
62
static void *hPROJMutex = NULL;
57
63
 
58
64
static projPJ       (*pfn_pj_init_plus)(const char *) = NULL;
59
65
static projPJ       (*pfn_pj_init)(int, char**) = NULL;
60
 
static projUV       (*pfn_pj_fwd)(projUV, projPJ) = NULL;
61
 
static projUV       (*pfn_pj_inv)(projUV, projPJ) = NULL;
62
66
static void     (*pfn_pj_free)(projPJ) = NULL;
63
67
static int      (*pfn_pj_transform)(projPJ, projPJ, long, int, 
64
68
                                    double *, double *, double * ) = NULL;
67
71
static char        *(*pfn_pj_get_def)(projPJ,int) = NULL;
68
72
static void         (*pfn_pj_dalloc)(void *) = NULL;
69
73
 
 
74
static projPJ (*pfn_pj_init_plus_ctx)( projCtx, const char * ) = NULL;
 
75
static int (*pfn_pj_ctx_get_errno)( projCtx ) = NULL;
 
76
static projCtx (*pfn_pj_ctx_alloc)(void) = NULL;
 
77
static void    (*pfn_pj_ctx_free)( projCtx ) = NULL;
 
78
 
70
79
#if (defined(WIN32) || defined(WIN32CE)) && !defined(__MINGW32__)
71
80
#  define LIBNAME      "proj.dll"
72
81
#elif defined(__MINGW32__)
112
121
    
113
122
    int         bCheckWithInvertProj;
114
123
    double      dfThreshold;
 
124
    
 
125
    projCtx     pjctx;
 
126
 
 
127
    int         InitializeNoLock( OGRSpatialReference *poSource, 
 
128
                                  OGRSpatialReference *poTarget );
 
129
 
 
130
    int         nMaxCount;
 
131
    double     *padfOriX;
 
132
    double     *padfOriY;
 
133
    double     *padfOriZ;
 
134
    double     *padfTargetX;
 
135
    double     *padfTargetY;
 
136
    double     *padfTargetZ;
115
137
 
116
138
public:
117
139
                OGRProj4CT();
165
187
#ifdef PROJ_STATIC
166
188
    pfn_pj_init = pj_init;
167
189
    pfn_pj_init_plus = pj_init_plus;
168
 
    pfn_pj_fwd = pj_fwd;
169
 
    pfn_pj_inv = pj_inv;
170
190
    pfn_pj_free = pj_free;
171
191
    pfn_pj_transform = pj_transform;
172
192
    pfn_pj_get_errno_ref = (int *(*)(void)) pj_get_errno_ref;
174
194
    pfn_pj_dalloc = pj_dalloc;
175
195
#if PJ_VERSION >= 446
176
196
    pfn_pj_get_def = pj_get_def;
177
 
#endif    
 
197
#endif
 
198
#if PJ_VERSION >= 480
 
199
    pfn_pj_ctx_alloc = pj_ctx_alloc;
 
200
    pfn_pj_ctx_free = pj_ctx_free;
 
201
    pfn_pj_init_plus_ctx = pj_init_plus_ctx;
 
202
    pfn_pj_ctx_get_errno = pj_ctx_get_errno;
 
203
#endif
178
204
#else
179
205
    CPLPushErrorHandler( CPLQuietErrorHandler );
180
206
 
187
213
 
188
214
    pfn_pj_init_plus = (projPJ (*)(const char *)) 
189
215
        CPLGetSymbol( pszLibName, "pj_init_plus" );
190
 
    pfn_pj_fwd = (projUV (*)(projUV,projPJ)) 
191
 
        CPLGetSymbol( pszLibName, "pj_fwd" );
192
 
    pfn_pj_inv = (projUV (*)(projUV,projPJ)) 
193
 
        CPLGetSymbol( pszLibName, "pj_inv" );
194
216
    pfn_pj_free = (void (*)(projPJ)) 
195
217
        CPLGetSymbol( pszLibName, "pj_free" );
196
218
    pfn_pj_transform = (int (*)(projPJ,projPJ,long,int,double*,
206
228
        CPLGetSymbol( pszLibName, "pj_get_def" );
207
229
    pfn_pj_dalloc = (void (*)(void*))
208
230
        CPLGetSymbol( pszLibName, "pj_dalloc" );
 
231
 
 
232
    /* PROJ 4.8.0 symbols */
 
233
    pfn_pj_ctx_alloc = (projCtx (*)( void ))
 
234
        CPLGetSymbol( pszLibName, "pj_ctx_alloc" );
 
235
    pfn_pj_ctx_free = (void (*)( projCtx ))
 
236
        CPLGetSymbol( pszLibName, "pj_ctx_free" );
 
237
    pfn_pj_init_plus_ctx = (projPJ (*)( projCtx, const char * ))
 
238
        CPLGetSymbol( pszLibName, "pj_init_plus_ctx" );
 
239
    pfn_pj_ctx_get_errno = (int (*)( projCtx ))
 
240
        CPLGetSymbol( pszLibName, "pj_ctx_get_errno" );
 
241
 
209
242
    CPLPopErrorHandler();
210
 
 
 
243
    CPLErrorReset();
211
244
#endif
212
245
 
 
246
    if (pfn_pj_ctx_alloc != NULL &&
 
247
        pfn_pj_ctx_free != NULL &&
 
248
        pfn_pj_init_plus_ctx != NULL &&
 
249
        pfn_pj_ctx_get_errno != NULL &&
 
250
        CSLTestBoolean(CPLGetConfigOption("USE_PROJ_480_FEATURES", "YES")))
 
251
    {
 
252
        CPLDebug("OGRCT", "PROJ >= 4.8.0 features enabled");
 
253
    }
 
254
    else
 
255
    {
 
256
        pfn_pj_ctx_alloc = NULL;
 
257
        pfn_pj_ctx_free = NULL;
 
258
        pfn_pj_init_plus_ctx = NULL;
 
259
        pfn_pj_ctx_get_errno = NULL;
 
260
    }
 
261
 
213
262
    if( pfn_pj_transform == NULL )
214
263
    {
215
264
        CPLError( CE_Failure, CPLE_AppDefined, 
237
286
{
238
287
    char        *pszNewProj4Def, *pszCopy;
239
288
    projPJ      psPJSource = NULL;
 
289
 
240
290
    CPLMutexHolderD( &hPROJMutex );
241
291
 
242
292
    if( !LoadProjLibrary() || pfn_pj_dalloc == NULL || pfn_pj_get_def == NULL )
333
383
{
334
384
    OGRProj4CT  *poCT;
335
385
 
336
 
    if( !LoadProjLibrary() )
 
386
    if( pfn_pj_init == NULL && !LoadProjLibrary() )
337
387
    {
338
388
        CPLError( CE_Failure, CPLE_NotSupported, 
339
389
                  "Unable to load PROJ.4 library (%s), creation of\n"
359
409
/*                   OCTNewCoordinateTransformation()                   */
360
410
/************************************************************************/
361
411
 
 
412
/**
 
413
 * Create transformation object.
 
414
 *
 
415
 * This is the same as the C++ function OGRCreateCoordinateTransformation().
 
416
 *
 
417
 * Input spatial reference system objects are assigned 
 
418
 * by copy (calling clone() method) and no ownership transfer occurs.
 
419
 *
 
420
 * OCTDestroyCoordinateTransformation() should
 
421
 * be used to destroy transformation objects. 
 
422
 *
 
423
 * The PROJ.4 library must be available at run-time.
 
424
 *
 
425
 * @param hSourceSRS source spatial reference system. 
 
426
 * @param hTargetSRS target spatial reference system. 
 
427
 * @return NULL on failure or a ready to use transformation object.
 
428
 */
 
429
 
362
430
OGRCoordinateTransformationH CPL_STDCALL 
363
431
OCTNewCoordinateTransformation(
364
432
    OGRSpatialReferenceH hSourceSRS, OGRSpatialReferenceH hTargetSRS )
386
454
    
387
455
    bCheckWithInvertProj = FALSE;
388
456
    dfThreshold = 0;
 
457
 
 
458
    nMaxCount = 0;
 
459
    padfOriX = NULL;
 
460
    padfOriY = NULL;
 
461
    padfOriZ = NULL;
 
462
    padfTargetX = NULL;
 
463
    padfTargetY = NULL;
 
464
    padfTargetZ = NULL;
 
465
 
 
466
    if (pfn_pj_ctx_alloc != NULL)
 
467
        pjctx = pfn_pj_ctx_alloc();
 
468
    else
 
469
        pjctx = NULL;
389
470
}
390
471
 
391
472
/************************************************************************/
407
488
            delete poSRSTarget;
408
489
    }
409
490
 
410
 
    CPLMutexHolderD( &hPROJMutex );
411
 
 
412
 
    if( psPJSource != NULL )
413
 
        pfn_pj_free( psPJSource );
414
 
 
415
 
    if( psPJTarget != NULL )
416
 
        pfn_pj_free( psPJTarget );
 
491
    if (pjctx != NULL)
 
492
    {
 
493
        pfn_pj_ctx_free(pjctx);
 
494
 
 
495
        if( psPJSource != NULL )
 
496
            pfn_pj_free( psPJSource );
 
497
 
 
498
        if( psPJTarget != NULL )
 
499
            pfn_pj_free( psPJTarget );
 
500
    }
 
501
    else
 
502
    {
 
503
        CPLMutexHolderD( &hPROJMutex );
 
504
 
 
505
        if( psPJSource != NULL )
 
506
            pfn_pj_free( psPJSource );
 
507
 
 
508
        if( psPJTarget != NULL )
 
509
            pfn_pj_free( psPJTarget );
 
510
    }
 
511
 
 
512
    CPLFree(padfOriX);
 
513
    CPLFree(padfOriY);
 
514
    CPLFree(padfOriZ);
 
515
    CPLFree(padfTargetX);
 
516
    CPLFree(padfTargetY);
 
517
    CPLFree(padfTargetZ);
417
518
}
418
519
 
419
520
/************************************************************************/
424
525
                            OGRSpatialReference * poTargetIn )
425
526
 
426
527
{
 
528
    if (pjctx != NULL)
 
529
    {
 
530
        return InitializeNoLock(poSourceIn, poTargetIn);
 
531
    }
 
532
 
427
533
    CPLMutexHolderD( &hPROJMutex );
428
 
 
 
534
    return InitializeNoLock(poSourceIn, poTargetIn);
 
535
}
 
536
 
 
537
/************************************************************************/
 
538
/*                         InitializeNoLock()                           */
 
539
/************************************************************************/
 
540
 
 
541
int OGRProj4CT::InitializeNoLock( OGRSpatialReference * poSourceIn, 
 
542
                                  OGRSpatialReference * poTargetIn )
 
543
 
 
544
{
429
545
    if( poSourceIn == NULL || poTargetIn == NULL )
430
546
        return FALSE;
431
547
 
539
655
        return FALSE;
540
656
    }
541
657
 
542
 
    psPJSource = pfn_pj_init_plus( pszProj4Defn );
 
658
    if (pjctx)
 
659
        psPJSource = pfn_pj_init_plus_ctx( pjctx, pszProj4Defn );
 
660
    else
 
661
        psPJSource = pfn_pj_init_plus( pszProj4Defn );
543
662
    
544
663
    if( psPJSource == NULL )
545
664
    {
546
 
        if( pfn_pj_get_errno_ref != NULL
 
665
        if( pjctx != NULL)
 
666
        {
 
667
            int pj_errno = pfn_pj_ctx_get_errno(pjctx);
 
668
 
 
669
            /* pfn_pj_strerrno not yet thread-safe in PROJ 4.8.0 */
 
670
            CPLMutexHolderD(&hPROJMutex);
 
671
            CPLError( CE_Failure, CPLE_NotSupported, 
 
672
                      "Failed to initialize PROJ.4 with `%s'.\n%s", 
 
673
                      pszProj4Defn, pfn_pj_strerrno(pj_errno) );
 
674
        }
 
675
        else if( pfn_pj_get_errno_ref != NULL
547
676
            && pfn_pj_strerrno != NULL )
548
677
        {
549
678
            int *p_pj_errno = pfn_pj_get_errno_ref();
588
717
        return FALSE;
589
718
    }
590
719
 
591
 
    psPJTarget = pfn_pj_init_plus( pszProj4Defn );
 
720
    if (pjctx)
 
721
        psPJTarget = pfn_pj_init_plus_ctx( pjctx, pszProj4Defn );
 
722
    else
 
723
        psPJTarget = pfn_pj_init_plus( pszProj4Defn );
592
724
    
593
725
    if( psPJTarget == NULL )
594
726
        CPLError( CE_Failure, CPLE_NotSupported, 
711
843
/* -------------------------------------------------------------------- */
712
844
/*      Do the transformation using PROJ.4.                             */
713
845
/* -------------------------------------------------------------------- */
714
 
    CPLMutexHolderD( &hPROJMutex );
 
846
    if (pjctx == NULL)
 
847
    {
 
848
        /* The mutex has already been created */
 
849
        CPLAssert(hPROJMutex != NULL);
 
850
        CPLAcquireMutex(hPROJMutex, 1000.0);
 
851
    }
715
852
        
716
853
    if (bCheckWithInvertProj)
717
854
    {
718
855
        /* For some projections, we cannot detect if we are trying to reproject */
719
856
        /* coordinates outside the validity area of the projection. So let's do */
720
857
        /* the reverse reprojection and compare with the source coordinates */
721
 
        
722
 
        double *ori_x = NULL;
723
 
        double *ori_y = NULL;
724
 
        double *ori_z = NULL;
725
 
        ori_x = (double*)CPLMalloc(sizeof(double)*nCount);
726
 
        memcpy(ori_x, x, sizeof(double)*nCount);
727
 
        ori_y = (double*)CPLMalloc(sizeof(double)*nCount);
728
 
        memcpy(ori_y, y, sizeof(double)*nCount);
 
858
        if (nCount > nMaxCount)
 
859
        {
 
860
            nMaxCount = nCount;
 
861
            padfOriX = (double*) CPLRealloc(padfOriX, sizeof(double)*nCount);
 
862
            padfOriY = (double*) CPLRealloc(padfOriY, sizeof(double)*nCount);
 
863
            padfOriZ = (double*) CPLRealloc(padfOriZ, sizeof(double)*nCount);
 
864
            padfTargetX = (double*) CPLRealloc(padfTargetX, sizeof(double)*nCount);
 
865
            padfTargetY = (double*) CPLRealloc(padfTargetY, sizeof(double)*nCount);
 
866
            padfTargetZ = (double*) CPLRealloc(padfTargetZ, sizeof(double)*nCount);
 
867
        }
 
868
        memcpy(padfOriX, x, sizeof(double)*nCount);
 
869
        memcpy(padfOriY, y, sizeof(double)*nCount);
729
870
        if (z)
730
871
        {
731
 
            ori_z = (double*)CPLMalloc(sizeof(double)*nCount);
732
 
            memcpy(ori_z, z, sizeof(double)*nCount);
 
872
            memcpy(padfOriZ, z, sizeof(double)*nCount);
733
873
        }
734
874
        err = pfn_pj_transform( psPJSource, psPJTarget, nCount, 1, x, y, z );
735
875
        if (err == 0)
736
876
        {
737
 
            double* target_x = (double*)CPLMalloc(sizeof(double)*nCount);
738
 
            double* target_y = (double*)CPLMalloc(sizeof(double)*nCount);
739
 
            double* target_z = NULL;
740
 
            memcpy(target_x, x, sizeof(double)*nCount);
741
 
            memcpy(target_y, y, sizeof(double)*nCount);
 
877
            memcpy(padfTargetX, x, sizeof(double)*nCount);
 
878
            memcpy(padfTargetY, y, sizeof(double)*nCount);
742
879
            if (z)
743
880
            {
744
 
                target_z = (double*)CPLMalloc(sizeof(double)*nCount);
745
 
                memcpy(target_z, z, sizeof(double)*nCount);
 
881
                memcpy(padfTargetZ, z, sizeof(double)*nCount);
746
882
            }
747
883
            
748
884
            err = pfn_pj_transform( psPJTarget, psPJSource , nCount, 1,
749
 
                                    target_x, target_y, target_z );
 
885
                                    padfTargetX, padfTargetY, (z) ? padfTargetZ : NULL);
750
886
            if (err == 0)
751
887
            {
752
888
                for( i = 0; i < nCount; i++ )
753
889
                {
754
890
                    if ( x[i] != HUGE_VAL && y[i] != HUGE_VAL &&
755
 
                        (fabs(target_x[i] - ori_x[i]) > dfThreshold ||
756
 
                         fabs(target_y[i] - ori_y[i]) > dfThreshold) )
 
891
                        (fabs(padfTargetX[i] - padfOriX[i]) > dfThreshold ||
 
892
                         fabs(padfTargetY[i] - padfOriY[i]) > dfThreshold) )
757
893
                    {
758
894
                        x[i] = HUGE_VAL;
759
895
                        y[i] = HUGE_VAL;
760
896
                    }
761
897
                }
762
898
            }
763
 
            
764
 
            CPLFree(target_x);
765
 
            CPLFree(target_y);
766
 
            CPLFree(target_z);
767
899
        }
768
 
        CPLFree(ori_x);
769
 
        CPLFree(ori_y);
770
 
        CPLFree(ori_z);
771
900
    }
772
901
    else
773
902
    {
787
916
 
788
917
        if( ++nErrorCount < 20 )
789
918
        {
 
919
            if (pjctx != NULL)
 
920
                /* pfn_pj_strerrno not yet thread-safe in PROJ 4.8.0 */
 
921
                CPLAcquireMutex(hPROJMutex, 1000.0);
 
922
 
790
923
            const char *pszError = NULL;
791
924
            if( pfn_pj_strerrno != NULL )
792
925
                pszError = pfn_pj_strerrno( err );
797
930
                          err );
798
931
            else
799
932
                CPLError( CE_Failure, CPLE_AppDefined, "%s", pszError );
 
933
 
 
934
            if (pjctx != NULL)
 
935
                /* pfn_pj_strerrno not yet thread-safe in PROJ 4.8.0 */
 
936
                CPLReleaseMutex(hPROJMutex);
800
937
        }
801
938
        else if( nErrorCount == 20 )
802
939
        {
805
942
                      err );
806
943
        }
807
944
 
 
945
        if (pjctx == NULL)
 
946
            CPLReleaseMutex(hPROJMutex);
808
947
        return FALSE;
809
948
    }
810
949
 
 
950
    if (pjctx == NULL)
 
951
        CPLReleaseMutex(hPROJMutex);
 
952
 
811
953
/* -------------------------------------------------------------------- */
812
954
/*      Potentially transform back to degrees.                          */
813
955
/* -------------------------------------------------------------------- */