2use,
INTRINSIC :: iso_c_binding
8TYPE(pj_object) :: pj, pjf, pjt
10TYPE(pj_coord_object) :: coordg, coordp, coordgc, &
11 vcoordg(4), vcoordp(4), vcoordgc(4)
13TYPE(pj_area_object) :: area
16CHARACTER(len=512) :: proj_stringf, proj_stringt
19proj_stringf =
'EPSG:4326'
20proj_stringt =
'+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs +type=crs'
21print*,
'from: ',trim(proj_stringf)
22print*,
'to: ',trim(proj_stringt)
25print*,
'Defining a projection object'
26pj = proj_create_crs_to_crs(pj_default_ctx, trim(proj_stringf)//char(0), &
27 trim(proj_stringt)//char(0), area)
29print*,
'Converting through a pj_coord object'
36coordp = proj_trans(pj, pj_fwd, coordg)
38coordgc = proj_trans(pj, pj_inv, coordp)
40print*,
'Original:',coordg
41print*,
'Projected:',coordp
42print*,
'Returned:',coordgc
45IF (abs(coordgc%x-coordg%x) > 1.0d-8 .OR. abs(coordgc%y-coordg%y) > 1.0d-8)
THEN
46 print*,
'Error pj_inv*pj_fwd failed or /= identity'
53pjf = proj_create(pj_default_ctx, trim(proj_stringf)//char(0))
54IF (.NOT.proj_associated(pjf))
THEN
55 print*,
'Failing to define projection object from ',trim(proj_stringf)
56 CALL print_error(pj_default_ctx)
59CALL print_info(proj_pj_info(pjf))
61pjt = proj_create(pj_default_ctx, trim(proj_stringt)//char(0))
62IF (.NOT.proj_associated(pjt))
THEN
63 print*,
'Failing to define projection object from ',trim(proj_stringt)
64 CALL print_error(pj_default_ctx)
67CALL print_info(proj_pj_info(pjt))
69print*,
'Defining a projection object'
70pj = proj_create_crs_to_crs_from_pj(pj_default_ctx, pjf, pjt, area, c_null_ptr)
71IF (.NOT.proj_associated(pjt))
THEN
72 print*,
'Failing to define projection object'
73 CALL print_error(pj_default_ctx)
77print*,
'Converting through a pj_coord array object'
80vcoordg(:)%x = (/45.0d0,45.0d0,46.0d0,46.0d0/)
81vcoordg(:)%y = (/11.0d0,12.0d0,11.0d0,12.0d0/)
85res = proj_trans_f(pj, pj_fwd, vcoordp)
87 print*,
'Failure in forward array transformation',res
88 CALL print_error(pj_default_ctx)
93res = proj_trans_f(pj, pj_inv, vcoordgc)
95 print*,
'Failure in inverse array transformation',res
96 CALL print_error(pj_default_ctx)
100print*,
'Original:',vcoordg%x,vcoordg%y
101print*,
'Projected:',vcoordp%x,vcoordp%y
102print*,
'Returned:',vcoordgc%x,vcoordgc%y
105IF (maxval(abs(vcoordgc%x-vcoordg%x)) > 1.0d-8 .OR. &
106 maxval(abs(vcoordgc%y-vcoordg%y)) > 1.0d-8)
THEN
107 print*,
'Error pj_inv*pj_fwd failed or /= identity'
115SUBROUTINE print_info(info)
116TYPE(pj_proj_info_object),
INTENT(in) :: info
118print*,
'==== Object information ===='
120print*,
'desc:',trim(
strtofchar(info%description,256))
121print*,
'def:',trim(
strtofchar(info%definition,256))
122print*,
'has inv:',info%has_inverse
123print*,
'accuracy:',info%accuracy
124print*,
'============================'
126END SUBROUTINE print_info
128SUBROUTINE print_error(ctx)
129TYPE(pj_context_object),
INTENT(in) :: ctx
131IF (proj_context_errno(ctx) /= 0)
THEN
133 print*,trim(
strtofchar(proj_errno_string(proj_context_errno(ctx)), 256))
136END SUBROUTINE print_error
138END PROGRAM proj_test6
Convert a null-terminated C string into a Fortran CHARACTER variable of the proper length.
Utility module for supporting Fortran 2003 C language interface module.
Fortran 2003 interface to the proj https://proj.org/ library, API version 6.