[cig-commits] [commit] devel, master: added the calculation of longitude and latitude of the observation grid for GMT display purposes when the observation grid is read back from a saved file on the disk (6fcc510)

[cig_noreply at geodynamics.org]

Repository : https://github.com/geodynamics/specfem3d_globe

On branches: devel,master
Link       : https://github.com/geodynamics/specfem3d_globe/compare/bc58e579b3b0838a0968725a076f5904845437ca...be63f20cbb6f462104e949894dbe205d2398cd7f

>---------------------------------------------------------------

commit 6fcc5103529f4d77a69b48f9976d4972d177ee9c
Author: Dimitri Komatitsch <komatitsch at lma.cnrs-mrs.fr>
Date:   Thu May 22 16:11:07 2014 +0200

    added the calculation of longitude and latitude of the observation grid for GMT display purposes when the observation grid is read back from a saved file on the disk


>---------------------------------------------------------------

6fcc5103529f4d77a69b48f9976d4972d177ee9c
 src/meshfem3D/compute_coordinates_grid.f90 | 23 ++++++++++++++++++++++-
 1 file changed, 22 insertions(+), 1 deletion(-)

diff --git a/src/meshfem3D/compute_coordinates_grid.f90 b/src/meshfem3D/compute_coordinates_grid.f90
index badbe03..e8b463c 100644
--- a/src/meshfem3D/compute_coordinates_grid.f90
+++ b/src/meshfem3D/compute_coordinates_grid.f90
@@ -378,7 +378,28 @@
     call bcast_all_dp(y_observation1D, NTOTAL_OBSERVATION)
     call bcast_all_dp(z_observation1D, NTOTAL_OBSERVATION)
 
-  else
+!   loop on all the chunks and then on all the observation nodes in each chunk
+    do ichunk = 1,NCHUNKS_MAX
+      do iy = 1,NY_OBSERVATION
+        do ix = 1,NX_OBSERVATION
+
+          x_top = x_observation(ix,iy,ichunk)
+          y_top = y_observation(ix,iy,ichunk)
+          z_top = z_observation(ix,iy,ichunk)
+
+          ! converts geocentric coordinates x/y/z to geographic radius/latitude/longitude (in degrees)
+          call xyz_2_rlatlon_dble(x_top,y_top,z_top,r,lat,lon)
+
+          ! store the values obtained for future display with GMT
+          if( lon > 180.0d0 ) lon = lon - 360.0d0
+          lon_observation(ix,iy,ichunk) = lon
+          lat_observation(ix,iy,ichunk) = lat
+
+        enddo
+      enddo
+    enddo
+
+  else ! of if(REUSE_EXISTING_OBSERVATION_SURF)
 
   ! for future GMT display
   if(myrank == 0) then