One approach is to change the code to add a patch object in this part of it (add ‘pv’ and ‘cv’ and the patch call):
set(polar(th,r1,'-k'),'LineWidth',1.2);set(polar(th,r2,'-k'),'LineWidth',1.2);
pv = linspace(0, 2*pi*1.01, 250); % Used in 'patch'
cv = [cos(pv); sin(pv)]; % Used in 'patch'
patch(cv(1,:), cv(2,:), 'w') % Plot 'patch'
title('$\eta_1=\eta_2=0.0$','Interpreter','latex','FontSize',9,'FontName','Times New Roman','FontWeight','Normal')
First, delete (or comment-out) these lines:
if ~isempty(I);
x(I,J) = 0;
y(I,J) = 0;
end
They create the centre square, and you need to get rid of it.
This plots everything, then creates a white patch object for the inside circle to cover the contour lines in it. (An alternative would be to get all the contour coordinates from ‘DH1’ and then delete the (x,y) pairs that are inside the inner circle, defined here by the patch coordinates ‘pv’ and ‘cv’. That could be done, however it would significantly slow the code. Creating the patch object is just easier.)
(The entire code is too long to post, exceeding the maximum character limit.)
EDIT — (2 Jan 2024 at 22:20)
The result is something like this —
.
