From 3f2867508e1a7c1fa55546b3471f1d84280139ce Mon Sep 17 00:00:00 2001 From: Luca De Siena Date: Tue, 4 Jun 2024 14:25:41 +0200 Subject: [PATCH] Changing to .png and commenting bin + Changed output format to png, tiff was occupying too much space + Cleaning descriptions in functions --- bin/Murat_body.m | 4 ++++ bin/Murat_codaCheck.m | 10 ++++++++-- bin/Murat_inversion.m | 4 ++-- bin/Murat_lsqlinQmean.m | 33 +++++++++++++++++++-------------- bin/Murat_plot.m | 14 +++++++------- bin/Murat_saveFigures.m | 8 ++++---- bin/Murat_saveFigures_2panels.m | 8 ++++---- 7 files changed, 48 insertions(+), 33 deletions(-) diff --git a/bin/Murat_body.m b/bin/Murat_body.m index 5734f3a..d3f4751 100644 --- a/bin/Murat_body.m +++ b/bin/Murat_body.m @@ -20,18 +20,22 @@ % energyRatioBodyCoda_i: energy ratio between body and coda waves % energyRatioCodaNoise_i: energy ratio between coda waves and noise +% Redefine windows with sampling for noise and direct waves int = bodyWindow*srate_i; lc = length(cursorCodaStart_i:cursorCodaEnd_i); cursor0 = floor(startNoise*srate_i); cursor0_1 = floor(cursor0 + int); cursor2 = floor(cursorPick_i + int-1); +% Direct wave normalised by window length spamp =... trapz(sp_i(cursorPick_i:cursor2,:))/int; +% Noise pre-pick normalised by window length spampn =... trapz(sp_i(cursor0:cursor0_1,:))/int; +% Coda normalised by window length spampc =... trapz(sp_i(cursorCodaStart_i:cursorCodaEnd_i,:))/lc; diff --git a/bin/Murat_codaCheck.m b/bin/Murat_codaCheck.m index 5a0df08..bc3778b 100644 --- a/bin/Murat_codaCheck.m +++ b/bin/Murat_codaCheck.m @@ -25,18 +25,22 @@ % cursorCodaStart_i: coda starting time after check on trace % cursorCodaEnd_i: coda end time after check on trace +%Define envelope duration t00 = tempis(1); lengthTempis = length(tempis); +% Method peak is generally valid for active seismics if isequal(peakDelayMethod,'Peak') tCoda_i =... (pktime_i-originTime_i)+peakDelay_i; - + +% Method constant is the most used in small local tomography elseif isequal(peakDelayMethod,'Constant') tCoda_i = tCm; - + +% Method travel is the standard at the regional scale elseif isequal(peakDelayMethod,'Travel') tCoda_i =... @@ -44,12 +48,14 @@ end +% Define the indexes along the seismogram cursorCodaStart_i =... floor((originTime_i - t00 + tCoda_i) * srate_i - 1); cursorCodaEnd_i =... floor(cursorCodaStart_i + tWm * srate_i - 1); +% Some seismograms are not long enough so you consider a shorter window if cursorCodaEnd_i > lengthTempis cursorCodaEnd_i = lengthTempis; end diff --git a/bin/Murat_inversion.m b/bin/Murat_inversion.m index db4e34a..9cb8ff7 100644 --- a/bin/Murat_inversion.m +++ b/bin/Murat_inversion.m @@ -130,7 +130,7 @@ modv_Qc(rcQc_k,4,k) = mQc; saveas(LcQc,fullfile(FPath, FLabel,'Tests/LCurve',FName)); - saveas(LcQc,fullfile(FPath, FLabel,'Tests/LCurve',FName),'tif'); + saveas(LcQc,fullfile(FPath, FLabel,'Tests/LCurve',FName),'png'); close(LcQc) %% @@ -181,7 +181,7 @@ modv_Q(rcQ_k,4,k) = mQ; saveas(LcCN,fullfile(FPath, FLabel,'Tests/LCurve',FName)); - saveas(LcCN,fullfile(FPath, FLabel,'Tests/LCurve',FName),'tif'); + saveas(LcCN,fullfile(FPath, FLabel,'Tests/LCurve',FName),'png'); close(LcCN) %% Checkerboards and spike inputs and checkerboard inversion diff --git a/bin/Murat_lsqlinQmean.m b/bin/Murat_lsqlinQmean.m index 434a758..931b9a2 100644 --- a/bin/Murat_lsqlinQmean.m +++ b/bin/Murat_lsqlinQmean.m @@ -3,30 +3,32 @@ % function [d1, const_Qc_k, constQmean_k, equationQ] =... % Murat_lsqlinQmean(tCm,tWm,Q_k,cf_k,sped,luntot_k,time0_k,rapsp_k) % -% INVERTS with weighted tikhonov and creates L-curve and data for Q. +% INVERTS with minimum least squares to obtain average Q % % Input parameters: -% tCm: lapse time -% tWm: coda window -% Q_k: inverse Qc +% tCm: starting coda time +% tWm: coda window length +% Q_k: average coda attenuation % cf_k: central frequency % sped: spectral decay -% luntot_k: total ray length -% time0_k: total travel time -% rapsp_k: energy ratio -% -% Output parameters: -% d1: data for average Q inversion -% const_Qc_k: constant for average Q inversion depending on Qc -% constQmean_k: average geometrical spreading and Q -% equationQ: sum of all three terms of the CN equation +% luntot_k: total length per frequency +% time0_k: travel time per frequency +% rapsp_k: energy ratio per frequency % +% Output parameters: +% d1: data for the inversion - variations from average +% const_Qc_k: constant obtained using the average source-station Qc +% constQmean_k: contains geometrical spreading, Q, uncertainties +% equationQ: equation to be compared with data in test + %% % Data creation for the true inversion, removing the pre-calculated % parameters. + +% Include info on Qc const_Qc_k = (tCm+tWm/2).^-sped... .*exp(-2*pi*Q_k*cf_k.*(tCm+tWm/2)); - +% Data vector for inversion of average Q d0 =... log(rapsp_k)/2/pi/cf_k + log(const_Qc_k)/2/pi/cf_k; @@ -35,13 +37,16 @@ cova = (G'*G)^(-1)*G'*cov(d0)*G*(G'*G)^(-1); +% Storing inverted parameters constQmean_k(:,1) = lsqlin(G,d0(:,1)); constQmean_k(:,2) = sqrt(diag(cova)); +% Calculate data vector for average Q and geometrical spreading d1 = d0 + ... constQmean_k(1,1)*log(luntot_k)/2/pi/cf_k... + time0_k*constQmean_k(2,1); +% Equation for average Q to fit data equationQ = -log(const_Qc_k)... -constQmean_k(1,1)*log(luntot_k)... -2*pi*cf_k*time0_k*constQmean_k(2,1); diff --git a/bin/Murat_plot.m b/bin/Murat_plot.m index 5fce938..58d5e54 100644 --- a/bin/Murat_plot.m +++ b/bin/Murat_plot.m @@ -71,7 +71,7 @@ storeFolder = 'Tests'; pathFolder =... fullfile(FPath,FLabel,storeFolder,FName_Cluster); - saveas(clustering,pathFolder,'tif'); + saveas(clustering,pathFolder,'png'); close(clustering) end @@ -105,7 +105,7 @@ ending,evestaz_pd,x,y,z,FName_peakDelay); pathFolder =... fullfile(FPath,FLabel,storeFolder,FName_peakDelay); - saveas(rays_peakDelay,pathFolder,'tif'); + saveas(rays_peakDelay,pathFolder,'png'); close(rays_peakDelay) %% @@ -117,7 +117,7 @@ Murat_imageRays(rma_Q,origin,ending,evestaz_Q,x,y,z,FName_Q); pathFolder =... fullfile(FPath, FLabel, storeFolder, FName_Q); - saveas(rays_Q,pathFolder,'tif'); + saveas(rays_Q,pathFolder,'png'); close(rays_Q) %% @@ -163,7 +163,7 @@ Qc_analysis = Murat_imageCheckQc(Qm_k,RZZ_k,... residualQc_k,luntot_Qc,Ac,sizeTitle,Qc_title,QcM); saveas(Qc_analysis, fullfile(FPath,FLabel,storeFolder,... - ['Qc_analysis_' fcName '_Hz']),'tif'); + ['Qc_analysis_' fcName '_Hz']),'png'); saveas(Qc_analysis, fullfile(FPath,FLabel,storeFolder,... ['Qc_analysis_' fcName '_Hz'])); close(Qc_analysis) @@ -179,7 +179,7 @@ pd_analysis = Murat_imageCheckPeakDelay(... time0PD,fitrobust_k,peakData_k,sizeTitle,pd_title); saveas(pd_analysis, fullfile(FPath,FLabel,storeFolder,... - ['PD_analysis_' fcName '_Hz']),'tif'); + ['PD_analysis_' fcName '_Hz']),'png'); saveas(pd_analysis, fullfile(FPath,FLabel,storeFolder,... ['PD_analysis_' fcName '_Hz'])); close(pd_analysis) @@ -208,7 +208,7 @@ Murat_imageCheckCN(equationQ,residualQ_k,d1,spreadAverageQ,... luntot_k,time0_k,energyRatio_k,A_k,Edirect_k,CN_title); saveas(CN_analysis, fullfile(FPath,FLabel,storeFolder,... - ['CN_analysis_' fcName '_Hz']),'tif'); + ['CN_analysis_' fcName '_Hz']),'png'); saveas(CN_analysis, fullfile(FPath,FLabel,storeFolder,... ['CN_analysis_' fcName '_Hz'])); close(CN_analysis) @@ -447,5 +447,5 @@ QcFrequency = Murat_imageQcFrequency(cf,... averageQcFrequency,sizeTitle,Qcf_title); FName = 'Qc_vs_frequency'; -saveas(QcFrequency, fullfile(FPath,FLabel,storeFolder,FName),'tif'); +saveas(QcFrequency, fullfile(FPath,FLabel,storeFolder,FName),'png'); close all diff --git a/bin/Murat_saveFigures.m b/bin/Murat_saveFigures.m index 1821eef..b272906 100644 --- a/bin/Murat_saveFigures.m +++ b/bin/Murat_saveFigures.m @@ -4,7 +4,7 @@ function Murat_saveFigures(figureName,Path) % % Input Parameters: % figureName: Matlab 3D plot -% Path: name of sections (tif) and figure (fig) +% Path: name of sections (png) and figure (fig) % % Output: % Three sections and 1 3D figure @@ -13,11 +13,11 @@ function Murat_saveFigures(figureName,Path) saveas(figureName,Path); view(90,0) ytickangle(45) -saveas(figureName,[Path '_SN'], 'tif'); +saveas(figureName,[Path '_SN'], 'png'); view(0,0) xtickangle(45) -saveas(figureName,[Path '_WE'], 'tif'); +saveas(figureName,[Path '_WE'], 'png'); view(0,90) -saveas(figureName,[Path '_H'], 'tif'); +saveas(figureName,[Path '_H'], 'png'); close(figureName) end \ No newline at end of file diff --git a/bin/Murat_saveFigures_2panels.m b/bin/Murat_saveFigures_2panels.m index 2190927..a03fba1 100644 --- a/bin/Murat_saveFigures_2panels.m +++ b/bin/Murat_saveFigures_2panels.m @@ -4,7 +4,7 @@ function Murat_saveFigures_2panels(figureName,Path) % % Input Parameters: % figureName: Matlab 3D plot -% Path: name of sections (tif) and figure (fig) +% Path: name of sections (png) and figure (fig) % % Output: % Three sections and 1 3D figure with 2 panels @@ -18,7 +18,7 @@ function Murat_saveFigures_2panels(figureName,Path) ax2 = subplot(1,2,2); view(ax2,90,0) ytickangle(45) -saveas(figureName,[Path '_SN'], 'tif'); +saveas(figureName,[Path '_SN'], 'png'); ax1 = subplot(1,2,1); view(ax1,0,0) @@ -26,13 +26,13 @@ function Murat_saveFigures_2panels(figureName,Path) ax2 = subplot(1,2,2); view(ax2,0,0) xtickangle(45) -saveas(figureName,[Path '_WE'], 'tif'); +saveas(figureName,[Path '_WE'], 'png'); ax1 = subplot(1,2,1); view(ax1,0,90) ax2 = subplot(1,2,2); view(ax2,0,90) -saveas(figureName,[Path '_H'], 'tif'); +saveas(figureName,[Path '_H'], 'png'); close(figureName)