- /
- 
        Ripples
        on 4 Nov 2024
        
        
 
    - 38
- 292
- 1
- 0
- 1059
 Cite your audio source here (if applicable): 
drawframe(1);drawframe(70);
 Write your drawframe function below
function drawframe(f)
persistent M W
if f==1
    % Capillary wave parameters
    s=.072; % water surface tension (N/m)
    % r=1025; % saltwater density (kg/m^3)
    r=1e3; % freshwater density (kg/m^3)
    a=.03; % ripple amplitude (m)
    T=1; % ripple period (s)
    w=2*pi/T; % angular frequency (rad/s)
    % Dispersion
    k=(w^2*r/s)^(1/3); % wavenumber (rad/m)
    L=2*pi/k; % wavelength (m)
    c=L/T; % wave celerity (m/s)
    ph=pi/2; % initial phase (rad)
    Fs=24; % Sampling frequency (96/4)
    t=1/Fs:1/Fs:4; % Time vector
    % Spatial domain
    x=-1:0.005:1;
    y=-1:0.005:1;
    % Water surface elevation calculation
    M=C(t,x,y,a,w,k,c,ph);
    % Water surface
    W=surf(x,y,squeeze(M(1,:,:)),'EdgeColor','none','FaceColor','b','FaceAlpha',1,'FaceLighting','gouraud');
    set(gca,'Position',[0 0 1 1])
    axis off
    axis([-.8 .8 -.8 .8 [-2 50]*max(M,[],'all')])    
    view(2)
    light(St='l',Po=[2 2 100*max(M,[],'all')],Col='c')
    lighting g
    material shiny
else
    % Updating the water surface
    W.ZData=squeeze(M(f,:,:));
end
end
function M=C(t,x,y,a,w,k,c,ph)
% Capillary wave calculation function
rng(0)
n=40; % Total number of ripples within the time duration
xc=-0.8+1.6*rand(1,n); yc=-0.8+1.6*rand(1,n); % Some random coordinates for each ripple
r=max(0.5,rand(1,length(xc))); % Some amplitude variance for each ripple
[X,Y]=meshgrid(x,y);
M=zeros(length(t),length(x),length(y),length(xc));
for j=1:numel(xc)
    for i=1:numel(t)
        tf=floor(length(t)/length(xc));
        tt=t(max(i-(j-1)*tf,1));
        if i>(j-1)*tf+1
            M(i,:,:,j)=a/length(xc)*(t(end-(j-1)*tf)-tt)/t(end-(j-1)*tf)*r(j)*sin(w*tt-k*sqrt((X-xc(j)).^2+(Y-yc(j)).^2)+ph);
        end
        x1=abs(x-xc(j))<c*tt*3;
        y1=abs(y-yc(j))<c*tt*3;
        % Mask
        m=x1.*y1';
        % 3D Hanning window mask to simulate ripple propagation & damping in space
        m(m==1)=ones(numel(find(x1==1)),numel(find(y1==1)))'.*(hann(numel(find(x1==1))).*hann(numel(find(y1==1)))')';
        M(i,:,:,j)=squeeze(M(i,:,:,j))'.*m';
    end
end
M=sum(M,4);
end


 

 
             
            