Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions tutorials/Tutorial_AlpineData.jl
Original file line number Diff line number Diff line change
Expand Up @@ -268,10 +268,10 @@ Vn = Vn*1000;
Vmagnitude = sqrt.(Ve.^2 + Vn.^2 + Vz.^2);

# ### 4.2 Interpolate topography on the grid
# At this stage we have the 3D velocity components on a grid. Yet, we don't have information yet about the elevation of the stations (as the provided data set did not give this).
# At this stage we have the 3D velocity components on a grid. Yet, we don't have information about the elevation of the stations (as the provided data set did not give this).
# We could ignore that and set the elevation to zero, which would allow saving the data directly.
# Yet, a better way is to load the topographic map of the area and interpolate the elevation to the velocity grid.
# As we have already the loaded the topographic map in section 1 of this tutorial, we can simply reuse it. To interpolate, we will use the function `interpolate_datafields_2D`
# As we already have loaded the topographic map in section 1 of this tutorial, we can simply reuse it. To interpolate, we will use the function `interpolate_datafields_2D`.

topo_v, fields_v = interpolate_datafields_2D(Topo, Lon, Lat)

Expand Down Expand Up @@ -319,7 +319,7 @@ Data_attribs = Dict(
"url"=>"https://doi.org/10.1029/2021JB023488",
"year"=>2022
)
# Now we are all set and can create a GeoData structure which along with metadata
# Now we are all set and can create a GeoData structure along with its metadata
Data = GeoData(Lon,Lat,Depth[:,:,:],(Vp=Vp[:,:,:],dVp=dlnVp[:,:,:]),Data_attribs);
# And then we save it again.
write_paraview(Data, "Rappisi2022")
Expand Down
6 changes: 3 additions & 3 deletions tutorials/Tutorial_Basic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ using GeophysicalModelGenerator;
Tomo_Alps_full = load_GMG("https://zenodo.org/records/10738510/files/Paffrath_2021_SE_Pwave.jld2?download=1")

# This is a so-called `GeoData` object, which is a 3D grid of seismic velocities as a function of longitude, latitude and depth, which can include various fields (here we only have a single field: `:dVp_Percentage`)
# We can save this in `VTK` format, which is a widely used format that can for exampke be read by the 3D open-source visualization tool [Paraview](https://www.paraview.org/):
# We can save this in `VTK` format, which is a widely used format that can for example be read by the 3D open-source visualization tool [Paraview](https://www.paraview.org/):
write_paraview(Tomo_Alps_full,"Tomo_Alps_full")

# We also uploaded a dataset with the topography of the Alpine region which can be downloaded with:
Expand All @@ -21,7 +21,7 @@ Topo_Alps = load_GMG("https://zenodo.org/records/10738510/files/AlpsTopo.jld2?do
# We can write this to disk as well
write_paraview(Topo_Alps,"Topo_Alps")

# If we open both datasets in Paraview, and changing both files from outline/solid colors to the corresponding data field, we see:
# If we open both datasets in Paraview, and change both files from outline/solid colors to the corresponding data field, we see:
# ![Basic_Tutorial_1](../assets/img/Basic_Tutorial_Paraview_1.png)
# Now we can change the colormap on the right side, marked by a red square. For topography we use the `Oleron` colormap, which you can download [here](https://www.fabiocrameri.ch/colourmaps/).
# For the tomography we use the `Roma` scientific colormap. You will now see a blue'ish box of the tomography, this is not the best color to visualise the data. Let's invert the colormap by clicking on the item marked by the blue arrow.
Expand Down Expand Up @@ -110,7 +110,7 @@ write_paraview(data_200km,"data_200km");

#
# ### 4. Cartesian data
# As you can see, the curvature or the Earth is taken into account here. Yet, for many applications it is more convenient to work in Cartesian coordinates (kilometers) rather then in geographic coordinates.
# As you can see, the curvature of the Earth is taken into account here. Yet, for many applications it is more convenient to work in Cartesian coordinates (kilometers) rather then in geographic coordinates.
# `GeophysicalModelGenerator` has a number of tools for this.
# First we need do define a `ProjectionPoint` around which we project the data
proj = ProjectionPoint(Lon=12.0,Lat =43)
Expand Down
2 changes: 1 addition & 1 deletion tutorials/Tutorial_GPS.jl
Original file line number Diff line number Diff line change
Expand Up @@ -140,7 +140,7 @@ height = interpol.(lon,lat)/1e3
# At this stage we have `lon/lat/height` of all points as well as velocity components
GPS_Sanchez_grid = GeoData(lon,lat,height,(Velocity_mm_year=(Ve,Vn,Vz),V_north=Vn*mm/yr, V_east=Ve*mm/yr, V_vertical=Vz*mm/yr, Vmagnitude = Vmagnitude*mm/yr, Topography = height*km))

# Save paraview is as always:
# Save to paraview as always:
write_paraview(GPS_Sanchez_grid, "GPSAlps_Sanchez_2017_grid")

# Opening and plotting the vertical field gives:
Expand Down
4 changes: 2 additions & 2 deletions tutorials/Tutorial_Jura.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@ using GeophysicalModelGenerator, GMT
Topo = import_topo(lat=[45.5,47.7], lon=[5, 8.1], file="@earth_relief_03s.grd")


# Next, we drape the geological map on top of the geological map.
# Next, we drape the geological map on top of the topography.
# The geological map was taken from the [2021 PhD thesis of Marc Schori](https://folia.unifr.ch/unifr/documents/313053) and saved as png map.
# We downloaded the pdf map, and cropped it to the lower left and upper right corners.
# The resulting map was uploaded to zenodo; it can be downloaded with
Expand All @@ -23,7 +23,7 @@ download_data("https://zenodo.org/records/10726801/files/SchoriM_Encl_01_Jura-ma
lowerleft = [4.54602510460251, 45.27456049638056, 0.0]
upperright = [8.948117154811715, 47.781282316442606, 0.0]

# We can now import the map with the `Screensho_To_GeoData` function:
# We can now import the map with the `screenshot_To_GeoData` function:
Geology = screenshot_to_GeoData("SchoriM_Encl_01_Jura-map_A1.png", lowerleft, upperright, fieldname=:geology_colors) # name should have "colors" in it

# You can "drape" this image on the topographic map with
Expand Down
2 changes: 1 addition & 1 deletion tutorials/Tutorial_LaPalma.jl
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@ using GeophysicalModelGenerator, GMT, DelimitedFiles
Topo = import_topo(lon = [-18.2, -17.5], lat=[28.4, 29.0], file="@earth_relief_15s.grd")

# Next, lets load the seismicity. The earthquake data is available on [https://www.ign.es/web/ign/portal/vlc-catalogo](https://www.ign.es/web/ign/portal/vlc-catalogo).
# We have filtered them and prepared a file with earthquake locations up to early November 2021 (from january 2021).
# We have filtered them and prepared a file with earthquake locations up to early November 2021 (from January 2021).
# Download that:
download_data("https://zenodo.org/records/10738510/files/EQ_events_all_info5_LaPalma_2021.dat","EQ_events_all_info5_LaPalma_2021.dat")
data_EQ = readdlm("EQ_events_all_info5_LaPalma_2021.dat")
Expand Down
4 changes: 2 additions & 2 deletions tutorials/Tutorial_NumericalModel_2D.jl
Original file line number Diff line number Diff line change
Expand Up @@ -129,7 +129,7 @@ Temp = fill(1350.0, nx, 1, nz);
add_box!(Phases, Temp, Grid2D; xlim=(-800,0.0), zlim=(-80.0, 0.0), phase = ConstantPhase(1));

# Next, we should define a `Trench` structure, which contains info about the trench which goes in 3D from `Start` - `End` coordinates (`x`,`y`)-coordinates respectively. As we are dealing with a 2D model, we set the `y`-coordinates to -100.0 and 100.0 respectively.
# Other parameters to be specified are `Thickness` (Slab thickness), `θ_max` (maximum slab dip angle), `Length` (length of slab), and `Lb` length of bending zoneof slab
# Other parameters to be specified are `Thickness` (Slab thickness), `θ_max` (maximum slab dip angle), `Length` (length of slab), and `Lb` (length of slab bending zone).

trench = Trench(Start=(0.0,-100.0), End=(0.0,100.0), Thickness=80.0, θ_max=45.0, Length=300, Lb=200, direction=-1.0);
add_slab!(Phases, Temp, Grid2D, trench, phase = ConstantPhase(1));
Expand Down Expand Up @@ -192,7 +192,7 @@ write_paraview(Grid2D,"Grid2D_SubductionCurvedOverriding");
# - `AddCylinder!`
#
# The help functions are quite self-explanatory, so we won't show it in detail here.
# If you have a topography surface or any other horizontal surface, you can surface with the cartesian grid with `above_surface` or `below_surface`.
# If you have a topography surface or any other horizontal surface, you can intersect the surface with the cartesian grid with `above_surface` or `below_surface`.
#
# Also, if you wish to take a seismic tomography as inspiration to set a slab geometry, you can interpolate it to a `CartGrid` with the same dimensions and use that with the julia `findall` function.

Expand Down
2 changes: 1 addition & 1 deletion tutorials/Tutorial_NumericalModel_3D.jl
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@ add_box!(Phases, Temp, Grid; xlim=(-800,0.0), ylim=(-400, 400.0), zlim=(-80.0, 0
T=SpreadingRateTemp(SpreadingVel=v_spreading_cm_yr, MORside="right"), StrikeAngle=30);


# And an an inclined part:
# And an inclined part:
AgeTrench = 800e3/(v_spreading_cm_yr/100)/1e6;
#add_box!(Phases, Temp, Grid; xlim=(0,300), ylim=(-400, 400.0), zlim=(-80.0, 0.0), phase = lith, Origin = (0,0,0), T=HalfspaceCoolingTemp(Age=0), DipAngle=0, StrikeAngle=30);
add_box!(Phases, Temp, Grid; xlim=(0,300), ylim=(-400, 400.0), zlim=(-80.0, 0.0), phase = lith, Origin = (0,0,0), T=HalfspaceCoolingTemp(Age=AgeTrench), DipAngle=30, StrikeAngle=30);
Expand Down
2 changes: 1 addition & 1 deletion tutorials/Tutorial_Votemaps.jl
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,7 @@ PSwave_Koulakov
# - Sum up the 3 models.
# - The resulting 3D map will have the number `3` at locations where all 3 models see a high-velocity anomaly (say a slab), implying that this feature is consistent between all models. Features that have a number 1 or 2 are only seen in a single or in 2/3 models.
#
# The result of this gives a feeling which features are consistent between the 3 models. It is ofcourse still subjective, as you still have to choose a criteria and as we are assuming in this that the 3 tomographic models are comparable (which they are not as they are produced using different amounts of datam, etc. etc.)
# The result of this gives a feeling which features are consistent between the 3 models. It is ofcourse still subjective, as you still have to choose a criteria and as we are assuming in this that the 3 tomographic models are comparable (which they are not as they are produced using different amounts of data, etc. etc.).
#
# So how do we create Votemaps?
# Doing this is rather simple:
Expand Down
Loading