Skip to content

Fix for geometric error in midpoint calculation and improvements in plot reproducibility #45

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
martinaavella wants to merge 1 commit into from
Closed

Fix for geometric error in midpoint calculation and improvements in plot reproducibility #45

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
4 changes: 2 additions & 2 deletions bin/Murat_checks.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -90,8 +90,8 @@

if availableVelocity == 0

qLat = mean(origin(1),ending(1));
qLon = mean(origin(2),ending(2));
qLat = mean([origin(1),ending(1)]);
qLon = mean([origin(2),ending(2)]);


gridPropagation.x = xM';
Expand Down
2 changes: 2 additions & 0 deletions bin/Murat_inversionQ.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -76,6 +76,8 @@
exitflag = 'Tikhonov';
output = [];
saveFigureAsImage(pathFolder);
savefig(gcf, [pathFolder '.fig']);
close(gcf);
fval = fval*obj0;

case 'Particle'
Expand Down
4 changes: 3 additions & 1 deletion bin/Murat_inversionQc.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,9 @@
Murat_tikhonovQc(PlotI,W*Ac_k,W*Qm_k,dampValue,x0.Qc);
eflag = 'Tikhonov';
output = [];
saveFigureAsImage(pathFolder);
saveFigureAsImage(pathFolder);
savefig(gcf, [pathFolder '.fig']);
close(gcf);
fval = fval*obj0;

case 'Particle'
Expand Down
43 changes: 28 additions & 15 deletions bin/Murat_plot.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -69,6 +69,21 @@
rtQck = retainQc(:,k);
rcQck = ray_crosses_Qc(:,k);

evst_pd = evst(rtpdk,:);
evst_Q = evst(rtQk,:);
Qm_k = Qm(rtQck,k);
RZZ_k = RZZ(rtQck,k);
avQcFreq(1,k) = sum(RZZ_k.*Qm_k)/sum(RZZ_k);
avQcFreq(2,k) = std(Qm_k);
modv_pd_k = modv_pd(:,:,k);
modv_Qc_k = modv_Qc(:,:,k);
modv_Q_k = modv_Q(:,:,k);
[X,Y,Z1,mPD]= Murat_fold(x,y,z,modv_pd_k(:,4));
[~,~,~,mQc] = Murat_fold(x,y,z,modv_Qc_k(:,4));
[~,~,~,mQ] = Murat_fold(x,y,z,modv_Q_k(:,4));
Z = Z1/1000;
evst_Qc = evst(rtQck,:);

if PlotRays == 1
% Murat_plot starts plotting the ray distribution if asked by the user.
% It stores the files in the corresponding folder.
Expand All @@ -77,21 +92,23 @@
% Peak Delay rays
FName_peakDelay = ['Rays_PeakDelay_' fcName '_Hz'];
rma_pd = rma(:,2:4,rtpdk)/1000;
evst_pd = evst(rtpdk,:);

rays_peakDelay = Murat_imageRays(rma_pd,origin,ending,evst_pd,...
x,y,z,FName_peakDelay);
saveFigureAsImage(makePath(storeFolder, FName_peakDelay));
close(rays_peakDelay)

drawnow limitrate

%%
% The next figure shows the rays for the total attenuation (Q)
FName_Q = ['Rays_Q_' fcName '_Hz'];
rma_Q = rma(:,2:4,rtQk)/1000;
evst_Q = evst(rtQk,:);

rays_Q = Murat_imageRays(rma_Q,origin,ending,evst_Q,x,y,z,...
FName_Q);
saveFigureAsImage(makePath(storeFolder, FName_Q));
close(rays_Q)
drawnow limitrate

end
% Tests
Expand All @@ -107,14 +124,12 @@

% Qc test
storeFolder = fullfile('Tests','Qc');
Qm_k = Qm(rtQck,k);
RZZ_k = RZZ(rtQck,k);

residualQc_k= residualQc(k);
luntot_Qc = luntot(rtQck)/1000;
Ac = Ac_i(rtQck,rcQck);

avQcFreq(1,k) = sum(RZZ_k.*Qm_k)/sum(RZZ_k);
avQcFreq(2,k) = std(Qm_k);


Qc_title = ['Qc check ' fcName ' Hz'];
Qc_analysis = Murat_imageCheckQc(Qm_k,RZZ_k,residualQc_k,...
Expand All @@ -123,6 +138,7 @@
saveFigureAsImage(p);
savefig(Qc_analysis, [p '.fig']);
close(Qc_analysis);
drawnow limitrate

% Peak delay test
storeFolder = fullfile('Tests','PeakDelay');
Expand All @@ -137,6 +153,7 @@
saveFigureAsImage(p);
savefig(pd_analysis, [p '.fig']);
close(pd_analysis);
drawnow limitrate

% Coda normalization test
storeFolder = fullfile('Tests','Q');
Expand All @@ -159,6 +176,7 @@
saveFigureAsImage(p);
savefig(CN_analysis, [p '.fig']);
close(CN_analysis);
drawnow limitrate

end

Expand All @@ -168,14 +186,7 @@
% work with the function "slice". All stored in the sub-folder.
storeFolder = fullfile('Results','PeakDelay');

modv_pd_k = modv_pd(:,:,k);
modv_Qc_k = modv_Qc(:,:,k);
modv_Q_k = modv_Q(:,:,k);
[X,Y,Z1,mPD]= Murat_fold(x,y,z,modv_pd_k(:,4));
[~,~,~,mQc] = Murat_fold(x,y,z,modv_Qc_k(:,4));
[~,~,~,mQ] = Murat_fold(x,y,z,modv_Q_k(:,4));
Z = Z1/1000;
evst_Qc = evst(rtQck,:);


% Peak delays results, using interpolation defined by 'divi'.
divi = 5;
Expand Down Expand Up @@ -350,6 +361,7 @@
modv_Qc_k,sTitle);
savefig(param_plot,makePath(storeFolder,FName_Param));
close(param_plot)
drawnow limitrate

% Use interpolated peakdelay and Qc
zi = (zi*1000)';
Expand All @@ -376,5 +388,6 @@
Murat_imageQcFrequency(cf, avQcFreq, sTitle, Qcf_title);
saveFigureAsImage(makePath('Results', 'Qc_vs_frequency'));
close all;
drawnow limitrate

end
95 changes: 91 additions & 4 deletions bin/saveFigureAsImage.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
% Resolve figure handle
if isempty(fig) || ~ishandle(fig)
fig = gcf;
fig.Color = 'white';
end

% Ensure filename ends with .png or treat as folder
Expand All @@ -37,15 +38,101 @@
end

% Use a fast exporter
origRenderer = fig.Renderer;
fig.Renderer = 'painter'; % good for raster PNGs
% origRenderer = fig.Renderer;
% fig.Renderer = 'painter'; % good for raster PNGs
ax = findall(fig, 'Type', 'axes');
for k = 1:numel(ax)
ax(k).Color = 'white'; % axes background
ax(k).XColor = 'k'; % tick and axis line color
ax(k).YColor = 'k';
if isprop(ax(k), 'ZColor')
ax(k).ZColor = 'k';
end
% Title and labels (if present)
if ~isempty(ax(k).Title)
ax(k).Title.Color = 'k';
end
if ~isempty(ax(k).XLabel)
ax(k).XLabel.Color = 'k';
end
if ~isempty(ax(k).YLabel)
ax(k).YLabel.Color = 'k';
end
end

% Make all text objects black (covers text(), annotation(), etc.)
txt = findall(fig, 'Type', 'text');
set(txt, 'Color', 'k');

% Configure all legends: white box background and black text
leg = findall(fig, 'Type', 'legend');
for k = 1:numel(leg)
% Legend text color
leg(k).TextColor = 'k';
% Legend background (opaque white)
if isprop(leg(k), 'Color')
leg(k).Color = 'white'; % older MATLAB also works
end
% For newer MATLAB, set the box face color to ensure opacity
if isprop(leg(k), 'BoxFace')
try
leg(k).BoxFace.ColorType = 'truecoloralpha'; % ensure supports alpha
leg(k).BoxFace.Color = [1 1 1]; % white
leg(k).BoxFace.FaceAlpha = 1; % fully opaque
catch
% ignore if property not supported
end
end
end

% Update all colorbars
cbs = findall(fig, 'Type', 'colorbar');
for k = 1:numel(cbs)
cb = cbs(k);
% Tick and label color (affects tick labels and tick marks)
if isprop(cb, 'Color')
cb.Color = 'k'; % black tick labels and tick marks
else
set(cb, 'Color', 'k');
end

% Label (colorbar title) – make sure it's black
if ~isempty(cb.Label)
cb.Label.Color = 'k';
end

% If legend-like box face exists (newer MATLAB), make it opaque white
if isprop(cb, 'BoxFace')
try
cb.BoxFace.ColorType = 'truecoloralpha';
cb.BoxFace.Color = [1 1 1]; % white
cb.BoxFace.FaceAlpha = 1; % fully opaque
catch
% ignore if not supported
end
else
% Fallback: try setting the colorbar background
try
cb.Color = 'k'; % keep ticks black
% no direct background property in older releases; use axes color
% set parent axes background white if needed:
ax = cb.Parent;
if isprop(ax, 'Color')
ax.Color = 'white';
end
catch
end
end
end


drawnow; % ensure content is up-to-date

% Use exportgraphics which is typically faster
% You can set resolution with 'Resolution',300 if needed.
exportgraphics(fig, outFile, 'BackgroundColor', 'white', 'Resolution', 72);
exportgraphics(fig, outFile, BackgroundColor='white', Resolution=72);

% restore renderer
fig.Renderer = origRenderer;
% fig.Renderer = origRenderer;

end
Loading