這是我寫的程式碼和圖形，但是就無法和原圖長的一樣，有人能告訴我少了什麼步驟嗎？？

r=0.0421;
q=0;
theta1=0.1; theta2=0.15;
sigma1=0.2; sigma2=0.25;
nu=1;
omega1=(1/nu)*( log(1-theta1*nu-sigma1*sigma1*nu/2));
m1=(r-q+omega1);　　%上面程式的裡面的參數m
omega2=(1/nu)*( log(1-theta2*nu-sigma2*sigma2*nu/2));
m2=(r-q+omega2);
%%%
c=zeros(99,99);
u=0.01:0.01:0.99;
for i=1:99
    F1=@(z1)quad(@(x)normcdf((z1-m1*1-theta1.*x)./(sigma1*sqrt(x))).*gampdf(x,1/nu,nu),0.001,100)-u(i);
    x1(i)=fzero(F1,u(i));
end
for j=1:99
    F2=@(z2)quad(@(x)normcdf((z2-m2*1-theta2.*x)./(sigma2*sqrt(x))).*gampdf(x,1/nu,nu),0.001,100)-u(j);
    x2(j)=fzero(F2,u(j));
end   
for i = 1:99
    for j = 1:99
        c(i,j)=quad(@(x)normcdf((x1(i)-m1*1-theta1.*x)./(sigma1*sqrt(x))).*normcdf((x2(j)-m2*1-theta2.*x)./(sigma2*sqrt(x))).*gampdf(x,1/nu,nu),0.001,100);
    end
end
contour(u,u,c)
axis([0, 1, 0, 1]);axis normal; %axis normal; axis square;
title('\theta_1=0.1,\theta_2=0.15,\sigma=0.2,\sigma=0.25,\nu=1');