How to plot this feather using a mathematical equation

I came across this image that intruiged me. It is a plot of a feather based on a mathematical equation. I attempted it, but couldn't get it right. Any idea how this could me done? I tried to use vectorisation method instead of a loop for now, but any method should be okay.

P.S: I'm still fairly new with MATLAB.

Here is the code I attempted

clear; close all; clc;
syms x y
f(x,y) = (-9/10) + ( 1 + (1/6) + ((1/(3*pi))*atan(30*( (3/50)-(y+(28/25)).^2 ))) )...
    .*(1-(((7*y)/10)+(2/10)).^10).*cos(x-(y/8)-(1/8)*(y+1).^2).^2000 ...
    +(1-sin((170*y)- 300*(x-(y/8)-(1/8)*(y+1).^2)^2 + (400/3)*(x-(y/8)-(1/8)*(y+1).^2)^3).^6)...
    .*((1/2)+((1/pi).*atan(400.*(1-(14.*(x-(1/20)-(y/8)-((1/8).*(y+1).^2).^2+(y.^2))).^2))));
figure()
fc = fcontour(f);
fc.LevelList = [0 0];

Would appreciate the help.

Answer 1

[Davide Masiello]

There:

clear,clc
[x,y]   = meshgrid(linspace(-0.4,0.9,1000),linspace(-1.5,1.1,1000));

xy = x-y/8-(1/8)*(y+1).^2;
A  = 1+1/6+atan(30*(3/50-(y+28/25).^2))/(3*pi);
B  = 1-(7*y/10+2/10).^10;
C  = 1-sin(170*y-300*xy.^2+(400/3)*xy.^3).^6;
D  = 1/2+atan(400*(1-(14*(xy-1/20).^2+y.^2).^2))/pi;

f = -9/10 + A.*B.*cos(xy).^2000 + C.*D;

f(f>0) = 1;
f(f<0) = 0;
contourf(x,y,f,[0 1])
axis equal
map = [1 1 1 ; 0 0 0];
colormap(map)

[Ahmed Mohamed Mansoor]

First of all thank you so much for the answer and quick response. I do have a couple of questions though.

1. How did you figure out the range for the meshgrid? Was it trial and error or was there a way?
2. Can you maybe explain what these lines actually do to the graph specifically

f(f>0) = 1;
f(f<0) = 0;

map = [1 1 1 ; 0 0 0];
colormap(map)

[Davide Masiello]

Hi @Ahmed Mohamed Mansoor.

1 - It was trial and error.

2 - From the picture you posted you can see that, in order to get the desired picture, you need only positive values of f.

Therefore, the lines

f(f>0) = 1;
f(f<0) = 0;

set any values of f higher than 0 equal to 1 and any value of f lower than 0 equal to 0.

This is useful because in the contour plot there are only two colours (black and white).

Accordingly, the lines

map = [1 1 1 ; 0 0 0];
colormap(map)

set a colormap with those two colours in it ([1 1 1] is white and [0 0 0] is black in rgb values).

[Ahmed Mohamed Mansoor]

Great!!. It makes much more sense now. I misunderstood the positive values of f(x,y) as only positive values of the range of x and y. But I understand it now and thank you again.

However, just to add on another request. Let's say i wanted to animate the line using the "drawnow" function on matlab. Would that be possible? If so how would you plot it out?

[Davide Masiello]

The code as it it now plots everything all at once.

How exactly do you want to animate it?

[Ahmed Mohamed Mansoor]

Yes you are correct.

Let's take a sin graph for example.

This is taken from the drawnow function documentation and it plots the points on a sing graph like an animation:

So would it possible to get it running like this?

clear; close all; clc;
h = animatedline;
axis([0 4*pi -1 1])
x = linspace(0,4*pi,2000);
for k = 1:length(x)
    y = sin(x(k));
    addpoints(h,x(k),y);
    drawnow
end