Added GMvD GMvMw scripts

etc/matlab/GroundMotionvsDistance.m

+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% Ground motion model (GMM) ground motion vs distance subplots
+%
+% Author: Peter Powers (pmpowers@usgs.gov)
+%         Demi Girot (dgirot@usgs.gov)
+%
+% Updated: 05/23/2023
+%
+% This script demonstrates how to call web services for the ground motion
+% models implemented in USGS National Sesmic Hazard Models and plot the
+% results. A JSON web service response is automatically converted to a
+% struct by the Matlab webread() function. This script shows how to plot
+% the GMMs requested and optionally the underlying models of epistemic
+% uncertainty.
+%
+% Users are STRONGLY encouraged to run the web services locally so as to
+% not overburden USGS servers. See the repository README for details on how
+% to do so: https://code.usgs.gov/ghsc/nshmp/nshmp-ws
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%% Plot a figure of GMM ground motion vs distance, sigmas, and epistemic 
+% uncertainty
+
+clearvars
+
+% The root url for the response spectra web service
+urlbase = "https://earthquake.usgs.gov/ws/nshmp/data/gmm/distance?";
+
+% If you are running nshmp-ws locally, use this URL:
+% urlbase = "http://localhost:8080/gmm/distance?";
+
+% Struct of ground motion model parameters
+in.Mw = 7.5;
+in.dip = 90;
+in.rake = 0;
+in.width = 14;
+in.rJB = 10;
+in.rMax = 300;
+in.rMin = 0.1;
+in.rRup = 10.3;
+in.rX = 10;
+in.vs30 = 760;
+in.vsInf = true;
+in.zHyp = 7.5;
+in.zTor = 0.5;
+
+% Uncomment to include any of these additional parameters in the struct
+% in.z2p5 = 1;
+% in.z1p0 = 1;
+% in.zSed = 1;
+
+figure
+plot_handles = [];
+
+% Set to true to show epistemic branches
+epi = false;
+
+gmms_2023 = ["ASK_14" "BSSA_14" "CB_14" "CY_14"];
+
+subplot(2,1,1)
+means = readGMvD(urlbase, gmms_2023, in).response.means.data;
+plotGmms(means, in, epi);
+
+subplot(2,1,2)
+sigmas = readGMvD(urlbase, gmms_2023, in).response.sigmas.data;
+plotSigmas(sigmas, in, epi);
+
+%% Shared Functions
+
+function plotGmms(gmm_data, in, epi)
+
+    plot_handles = [];
+    legend_labels = {};
+
+    % loop means response array
+    for j = 1 : length(gmm_data) 
+    
+        gmm_xs = gmm_data(j).data.xs;
+        gmm_ys = gmm_data(j).data.ys;
+        gmm_tree = gmm_data(j).tree;
+
+        % if 0.01 absent use PGA (0.001 is used for PGA)
+        if gmm_xs(1) == 0.001 && gmm_xs(2) ~= 0.01
+            gmm_xs(1) = 0.01;
+        end
+
+        % Plot the total spectrum
+        ph = loglog(gmm_xs, gmm_ys, ...
+            'LineWidth',3, ...
+            'LineStyle','-');
+        
+        color = get(ph,'Color');
+        hold on;
+        
+        plot_handles = [plot_handles ph];
+        legend_labels{end+1} = gmm_data(j).label;
+        
+        % Plot epistemic spectra, if present
+        if epi && ~isempty(gmm_tree)
+            for k = 1 : length(gmm_tree)
+                id = gmm_tree(k).id;
+                epi_ys = gmm_tree(k).values;
+                
+                phe = loglog(gmm_xs, epi_ys, ...
+                    'LineWidth',1, ...
+                    'LineStyle',':', ...
+                    'Color',color); 
+
+                plot_handles = [plot_handles phe];
+                legend_labels{end+1} = [gmm_data(j).label ' ' id];
+            end
+        end
+    end
+    
+    grid on
+
+    xlabel('Distance (km)','FontSize',16)
+    ylabel('Median Ground Motion (g)','FontSize',16)
+    title(['M',num2str(in.Mw),', ',num2str(in.rRup),'km'],'FontSize',20)
+    xlim([0.1 300]);
+
+    set(gca,'FontSize',20);
+    set(gca,'XTick',[0.1 1 10 100]);
+    set(gca,'XTickLabel',{'0.1','1','10','100'});
+    
+    legend(plot_handles,legend_labels,'Location','best');
+end
+
+function plotSigmas(gmm_data, in, epi)
+
+    plot_handles = [];
+    legend_labels = {};
+
+    % loop sigmas response array
+    for m = 1 : length(gmm_data) 
+    
+        sig_xs = gmm_data(m).data.xs;
+        sig_ys = gmm_data(m).data.ys;
+        sig_tree = gmm_data(m).tree;
+
+        % if 0.01 absent use PGA (0.001 is used for PGA)
+        if sig_xs(1) == 0.001 && sig_xs(2) ~= 0.01
+            sig_xs(1) = 0.01;
+        end
+
+        % Plot the total spectrum
+        ph = loglog(sig_xs, sig_ys, ...
+            'LineWidth',3, ...
+            'LineStyle','-');
+        
+        color = get(ph,'Color');
+        hold on;
+        
+        plot_handles = [plot_handles ph];
+        legend_labels{end+1} = gmm_data(m).label;
+        
+        % Plot epistemic spectra, if present
+        if epi && ~isempty(sig_tree)
+            for n = 1 : length(sig_tree)
+                id = sig_tree(n).id;
+                epi_ys = sig_tree(n).values;
+                
+                phe = loglog(sig_xs, epi_ys, ...
+                    'LineWidth',1, ...
+                    'LineStyle',':', ...
+                    'Color',color); 
+
+                plot_handles = [plot_handles phe];
+                legend_labels{end+1} = [gmm_data(m).label ' ' id];
+            end
+        end
+    end
+    
+    grid on
+
+    xlabel('Distance (km)','FontSize',16)
+    ylabel('Standard Deviation','FontSize',16)
+    title(['M',num2str(in.Mw),', ',num2str(in.rRup),'km'],'FontSize',20)
+    xlim([0.1 300]);
+
+    set(gca,'FontSize',20);
+    set(gca,'XTick',[0.1 1 10 100]);
+    set(gca,'XTickLabel',{'0.1','1','10','100'})
+    
+    legend(plot_handles,legend_labels,'Location','best');
+end
+
+% Read reponse spectra from a NSHM web service for the supplied
+% ground motion models (GMMs) and GMM input parameters
+function response = readGMvD(urlbase, gmms, in)
+
+url = urlbase;
+for i = 1 : size(gmms, 2)
+    if i == 1
+        url = url + "gmm=" + gmms(i);
+    else
+        url = url + "&gmm=" + gmms(i);
+    end
+end
+url = url + ...
+    "&Mw=" + num2str(in.Mw) + ...
+    "&dip=" + num2str(in.dip) + ...
+    "&rake=" + num2str(in.rake) + ...
+    "&width=" + num2str(in.width) + ...
+    "&rJB=" + num2str(in.rJB) + ...
+    "&rMax=" + num2str(in.rMax) + ...
+    "&rMin=" + num2str(in.rMin) + ...
+    "&rRup=" + num2str(in.rRup) + ...
+    "&rX=" + num2str(in.rX) + ...
+    "&vs30=" + num2str(in.vs30) + ...
+    "&vsInf=" + string(in.vsInf) + ...
+    "&zHyp=" + num2str(in.zHyp) + ...
+    "&zTop=" + num2str(in.zTor);
+
+if isfield(in, 'z2p5')
+    url = url + num2str(in.z2p5);
+end
+
+if isfield(in, 'z1p0')
+    url = url + num2str(in.z1p0);
+end
+
+if isfield(in, 'zSed')
+    url = url + num2str(in.zSed);
+end
+
+response = webread(url);
+
+end
\ No newline at end of file

etc/matlab/GroundMotionvsMagnitude.m

+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+% Ground motion model (GMM) ground motion vs magnitude subplots
+%
+% Author: Peter Powers (pmpowers@usgs.gov)
+%         Demi Girot (dgirot@usgs.gov)
+%
+% Updated: 05/23/2023
+%
+% This script demonstrates how to call web services for the ground motion
+% models implemented in USGS National Sesmic Hazard Models and plot the
+% results. A JSON web service response is automatically converted to a
+% struct by the Matlab webread() function. This script shows how to plot
+% the GMMs requested and optionally the underlying models of epistemic
+% uncertainty.
+%
+% Users are STRONGLY encouraged to run the web services locally so as to
+% not overburden USGS servers. See the repository README for details on how
+% to do so: https://code.usgs.gov/ghsc/nshmp/nshmp-ws
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%% Plot a figure of GMM ground motion vs magnitude, sigmas, and epistemic 
+% uncertainty
+
+clearvars
+
+% The root url for the response spectra web service
+urlbase = "https://earthquake.usgs.gov/ws/nshmp/data/gmm/magnitude?";
+
+% If you are running nshmp-ws locally, use this URL:
+% urlbase = "http://localhost:8080/gmm/magnitude?";
+
+% Struct of ground motion model parameters
+in.Mw = 7.5;
+in.dip = 90;
+in.rake = 0;
+in.width = 14;
+in.rJB = 10;
+in.rMax = 300;
+in.rMin = 0.1;
+in.rRup = 10.3;
+in.rX = 10;
+in.vs30 = 760;
+in.vsInf = true;
+in.zHyp = 7.5;
+in.zTor = 0.5;
+
+% Uncomment to include any of these additional parameters in the struct
+% in.z2p5 = 1;
+% in.z1p0 = 1;
+% in.zSed = 1;
+
+figure
+plot_handles = [];
+
+% Set to true to show epistemic branches
+epi = false;
+
+gmms_2023 = ["ASK_14" "BSSA_14" "CB_14" "CY_14"];
+
+subplot(2,1,1)
+means = readGMvMw(urlbase, gmms_2023, in).response.means.data;
+plotGmms(means, in, epi);
+
+subplot(2,1,2)
+sigmas = readGMvMw(urlbase, gmms_2023, in).response.sigmas.data;
+plotSigmas(sigmas, in, epi);
+
+%% Shared Functions
+
+function plotGmms(gmm_data, in, epi)
+
+    plot_handles = [];
+    legend_labels = {};
+
+    % loop means response array
+    for j = 1 : length(gmm_data) 
+    
+        gmm_xs = gmm_data(j).data.xs;
+        gmm_ys = gmm_data(j).data.ys;
+        gmm_tree = gmm_data(j).tree;
+
+        % if 0.01 absent use PGA (0.001 is used for PGA)
+        if gmm_xs(1) == 0.001 && gmm_xs(2) ~= 0.01
+            gmm_xs(1) = 0.01;
+        end
+
+        % Plot the total spectrum
+        ph = loglog(gmm_xs, gmm_ys, ...
+            'LineWidth',3, ...
+            'LineStyle','-');
+        
+        color = get(ph,'Color');
+        hold on;
+        
+        plot_handles = [plot_handles ph];
+        legend_labels{end+1} = gmm_data(j).label;
+        
+        % Plot epistemic spectra, if present
+        if epi && ~isempty(gmm_tree)
+            for k = 1 : length(gmm_tree)
+                id = gmm_tree(k).id;
+                epi_ys = gmm_tree(k).values;
+                
+                phe = loglog(gmm_xs, epi_ys, ...
+                    'LineWidth',1, ...
+                    'LineStyle',':', ...
+                    'Color',color); 
+
+                plot_handles = [plot_handles phe];
+                legend_labels{end+1} = [gmm_data(j).label ' ' id];
+            end
+        end
+    end
+    
+    grid on
+
+    xlabel('Distance (km)','FontSize',16)
+    ylabel('Median Ground Motion (g)','FontSize',16)
+    title(['M',num2str(in.Mw),', ',num2str(in.rRup),'km'],'FontSize',20)
+    set(gca,'FontSize',20);
+    
+    legend(plot_handles,legend_labels,'Location','best');
+end
+
+function plotSigmas(gmm_data, in, epi)
+
+    plot_handles = [];
+    legend_labels = {};
+
+    % loop sigmas response array
+    for m = 1 : length(gmm_data) 
+    
+        sig_xs = gmm_data(m).data.xs;
+        sig_ys = gmm_data(m).data.ys;
+        sig_tree = gmm_data(m).tree;
+
+        % if 0.01 absent use PGA (0.001 is used for PGA)
+        if sig_xs(1) == 0.001 && sig_xs(2) ~= 0.01
+            sig_xs(1) = 0.01;
+        end
+
+        % Plot the total spectrum
+        ph = loglog(sig_xs, sig_ys, ...
+            'LineWidth',3, ...
+            'LineStyle','-');
+        
+        color = get(ph,'Color');
+        hold on;
+        
+        plot_handles = [plot_handles ph];
+        legend_labels{end+1} = gmm_data(m).label;
+        
+        % Plot epistemic spectra, if present
+        if epi && ~isempty(sig_tree)
+            for n = 1 : length(sig_tree)
+                id = sig_tree(n).id;
+                epi_ys = sig_tree(n).values;
+                
+                phe = loglog(sig_xs, epi_ys, ...
+                    'LineWidth',1, ...
+                    'LineStyle',':', ...
+                    'Color',color); 
+
+                plot_handles = [plot_handles phe];
+                legend_labels{end+1} = [gmm_data(m).label ' ' id];
+            end
+        end
+    end
+    
+    grid on
+
+    xlabel('Distance (km)','FontSize',16)
+    ylabel('Standard Deviation','FontSize',16)
+    title(['M',num2str(in.Mw),', ',num2str(in.rRup),'km'],'FontSize',20)
+    set(gca,'FontSize',20);
+    
+    legend(plot_handles,legend_labels,'Location','best');
+end
+
+% Read reponse spectra from a NSHM web service for the supplied
+% ground motion models (GMMs) and GMM input parameters
+function response = readGMvMw(urlbase, gmms, in)
+
+url = urlbase;
+for i = 1 : size(gmms, 2)
+    if i == 1
+        url = url + "gmm=" + gmms(i);
+    else
+        url = url + "&gmm=" + gmms(i);
+    end
+end
+url = url + ...
+    "&Mw=" + num2str(in.Mw) + ...
+    "&dip=" + num2str(in.dip) + ...
+    "&rake=" + num2str(in.rake) + ...
+    "&width=" + num2str(in.width) + ...
+    "&rJB=" + num2str(in.rJB) + ...
+    "&rMax=" + num2str(in.rMax) + ...
+    "&rMin=" + num2str(in.rMin) + ...
+    "&rRup=" + num2str(in.rRup) + ...
+    "&rX=" + num2str(in.rX) + ...
+    "&vs30=" + num2str(in.vs30) + ...
+    "&vsInf=" + string(in.vsInf) + ...
+    "&zHyp=" + num2str(in.zHyp) + ...
+    "&zTop=" + num2str(in.zTor);
+
+if isfield(in, 'z2p5')
+    url = url + num2str(in.z2p5);
+end
+
+if isfield(in, 'z1p0')
+    url = url + num2str(in.z1p0);
+end
+
+if isfield(in, 'zSed')
+    url = url + num2str(in.zSed);
+end
+
+response = webread(url);
+
+end
\ No newline at end of file