From c328123c0461695a42ff9997050a784a1f93b3dd Mon Sep 17 00:00:00 2001 From: Jacob Frasunkiewicz Date: Thu, 5 Jun 2025 10:41:47 +0200 Subject: [PATCH 1/5] new function to import tomographies from NetCDF with tests --- src/data_import.jl | 30 +++++++++++++++++++++++++++++- test/test_data_import.jl | 13 +++++++++++++ 2 files changed, 42 insertions(+), 1 deletion(-) diff --git a/src/data_import.jl b/src/data_import.jl index 72952caa..72a4f404 100644 --- a/src/data_import.jl +++ b/src/data_import.jl @@ -6,7 +6,7 @@ using LightXML -export screenshot_to_GeoData, screenshot_to_CartData, screenshot_to_UTMData, getlonlatdepthmag_QuakeML +export screenshot_to_GeoData, screenshot_to_CartData, screenshot_to_UTMData, getlonlatdepthmag_QuakeML, tomo_2_GeoData # import CSV data using standard library functions # here we assume that the data is indeed comma separated and that comments are preceded with a "#" @@ -350,3 +350,31 @@ function getlonlatdepthmag_QuakeML(filename::String) Data_ISC = GeoData(lon, lat, -1 * depth / 1.0e3, (Magnitude = mag, Depth = -1 * depth / 1.0e3 * km)) return Data_ISC end + +""" + Read_TomoData(filename::String) + +Reads a seismic tomography dataset from a NetCDF file as a GeoData object. The keyword argument `vel_type::String` allows you to specify the type of velocity data to extract (default is "vs" for shear wave velocity). +The function assumes that the NetCDF file contains variables for depth, longitude, latitude, and the specified velocity type. + +tomodata = tomo_2_GeoData("path/to/tomo_data.nc") +""" +function tomo_2_GeoData(filename::String; vel_type::String = "vs") + + # Open the NetCDF file + data = NCDataset(filename) + + # Extract the variables + depth = data["depth"][:] + lon = data["longitude"][:] + lat = data["latitude"][:] + vel = data[vel_type][:,:,:] + + # create lon, lat, depth grid + Lon, Lat, Depth = lonlatdepth_grid(lon, lat, .- depth) + + # greate GeoDate struct + Tomo_data = GeoData(Lon, Lat, Depth, (vel=vel,)) + + return Tomo_data +end \ No newline at end of file diff --git a/test/test_data_import.jl b/test/test_data_import.jl index 07f27f76..3b716068 100644 --- a/test/test_data_import.jl +++ b/test/test_data_import.jl @@ -87,3 +87,16 @@ data_Image = screenshot_to_UTMData(filename, Corner_LowerLeft, Corner_UpperRight Data_ISC = getlonlatdepthmag_QuakeML("test_files/ISCTest.xml"); @test Value(Data_ISC.depth[1]) == -13.0km @test Data_ISC.fields.Magnitude[1] == 5.8 + +# test the import of a seismic tomography model from a NetCDF file +tomo_file = "test_files/test_tomo_data.nc" +Tomo_data = tomo_2_GeoData(tomo_file) +test_array = [i for i in 0:1:10] + +@test Tomo_data.lon.val[:,1,1] ≈ test_array +@test Tomo_data.lat.val[1,:,1] ≈ test_array +@test Tomo_data.depth.val[1,1,:] ≈ .- test_array + +for i in eachindex(test_array) + @test Tomo_data.fields.vel[i,i,i] ≈ i +end \ No newline at end of file From fa53244eb1bb8904ebb3a7aedaa01c7c3a2ed4fa Mon Sep 17 00:00:00 2001 From: Jacob Frasunkiewicz <94373204+jfrasunk@users.noreply.github.com> Date: Thu, 5 Jun 2025 10:44:50 +0200 Subject: [PATCH 2/5] Add files via upload small test file --- test/test_files/test_tomo_data.nc | Bin 0 -> 18409 bytes 1 file changed, 0 insertions(+), 0 deletions(-) create mode 100644 test/test_files/test_tomo_data.nc diff --git a/test/test_files/test_tomo_data.nc b/test/test_files/test_tomo_data.nc new file mode 100644 index 0000000000000000000000000000000000000000..5fa0aaeb491c19769947610dd531cacb1254a92f GIT binary patch literal 18409 zcmeI3PiP}`6vuy)Nn^G<-B!1Yt9EsaRZy4E?nc2w*@QHkMB8+2TJRzzbf)bc!h0iZ>67c(CY24<1wmMGyopDyw)AJ?zPw;G$Pg>+d`7_u5$7LuA?R!27|>H}m`R z=FR8#es5-ycBeLdd1&y#!BjS@xZE=3_bSy)n1_oV zDP^HcW(<2PjXlbV^{~AfM6Y`OpH!hnM)2WR)8iFWduB|b#!jh&miajI5V~4v+lwvb~NUw z2eHC{^E8famC?@XzU7%xv#Cz&2EfGiPuFKFjb?qWvD7S=W-7`VN4>2EXh5ShTZ#J* zn*Mv!5}USl&r{k%ZJIrt4mJYs-@fUgR37Ex=RZI1eoK7SZRT4apc&jIQFkPOfunD}F+@G^^am3^9Ttf&y8Dr8&6R5}I2BwD zLO5O>OO5iCPPplJqChV%tdgAGFks8?bc3)xb>VzFUWZeK{B!vW=hwWIXQn3elZE_6 zTIbhIx-021alv4fzxZ(P71<`KuiwVOXo^mM~AekwVF#yyW@asdu##bI!S6<{nLc+V?6@7jr5_Gf-L4;WH-cq13OB;lPHS`Rb$2_6*4*lB_0jjo zLp@ZA5D)@)lt8sSHxp+iBvoBF>zgN6#kU~7F>yV=rf{C%dVVh=Yf61M%%4PsfDjM@ zLO=)z0U;m+gn$qb0zyCt2!XqUz+ab)A4Srv(GjDgM$Z}@GwK@6867vOOi7Ae_EVH5 zPB;h$5J*}as3OgqsBnl0NIXy=B#Z+;e~8lb2?yZ-0!fPlet!_9_!SPq0R)m32mJmo zO7SZkgaZg9Ee`nmJ5h>X;UFA9AZc;H-@l7e{0ax*00K#i1OEL?l;T%72nP^IS{&qx z@LKsrm DZaLsd literal 0 HcmV?d00001 From aff0823510f73ca9e6ba4066c968584bd743212f Mon Sep 17 00:00:00 2001 From: Jacob Frasunkiewicz Date: Thu, 5 Jun 2025 15:32:05 +0200 Subject: [PATCH 3/5] adapted test of import_topo to region with more stable topography (Mainz) --- test/test_GMT.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/test/test_GMT.jl b/test/test_GMT.jl index 612c7337..9528ef3d 100644 --- a/test/test_GMT.jl +++ b/test/test_GMT.jl @@ -1,11 +1,11 @@ using Test using GeophysicalModelGenerator, GMT -Topo = import_topo(lat = [30, 31], lon = [50, 51]) -@test sum(Topo.depth.val) ≈ 2777.5705 +Topo = import_topo(lon = [8, 9], lat = [50, 51]) +@test sum(Topo.depth.val) ≈ 1076.7045 rtol = 1e-1 -Topo = import_topo([50, 51, 30, 31]); -@test sum(Topo.depth.val) ≈ 2777.5705 +Topo = import_topo([8, 9, 50, 51]); +@test sum(Topo.depth.val) ≈ 1076.7045 rtol = 1e-1 test_fwd = import_GeoTIFF("test_files/length_fwd.tif", fieldname = :forward) @test maximum(test_fwd.fields.forward) ≈ 33.17775km From 62c478dab3f8cb5fa1e2c1f051b81dd9f3166c17 Mon Sep 17 00:00:00 2001 From: Jacob Frasunkiewicz Date: Thu, 5 Jun 2025 15:35:30 +0200 Subject: [PATCH 4/5] added documentation --- docs/src/man/dataimport.md | 1 + 1 file changed, 1 insertion(+) diff --git a/docs/src/man/dataimport.md b/docs/src/man/dataimport.md index fe6c95ff..d4d7098d 100644 --- a/docs/src/man/dataimport.md +++ b/docs/src/man/dataimport.md @@ -8,4 +8,5 @@ GeophysicalModelGenerator.screenshot_to_CartData GeophysicalModelGenerator.screenshot_to_UTMData GeophysicalModelGenerator.import_topo GeophysicalModelGenerator.import_GeoTIFF +GeophysicalModelGenerator.tomo_2_GeoData ``` From 2fb72304ea3d76d377774cebee200eda5350530f Mon Sep 17 00:00:00 2001 From: Jacob Frasunkiewicz Date: Thu, 5 Jun 2025 15:36:54 +0200 Subject: [PATCH 5/5] fixed typos --- src/data_import.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/data_import.jl b/src/data_import.jl index 72a4f404..e983f305 100644 --- a/src/data_import.jl +++ b/src/data_import.jl @@ -373,7 +373,7 @@ function tomo_2_GeoData(filename::String; vel_type::String = "vs") # create lon, lat, depth grid Lon, Lat, Depth = lonlatdepth_grid(lon, lat, .- depth) - # greate GeoDate struct + # create GeoData struct Tomo_data = GeoData(Lon, Lat, Depth, (vel=vel,)) return Tomo_data