等高線を繋いでみる

等高線を描いたら，なんか領域の境目でブチギレちまってだせえ．そんな時には，自分で繋いでしまいましょう．

#サンプルデータ
X=linspace(0,1,11); #[0,1]を10糖分
Y=linspace(-1,1,21);#[-1,1]を20糖分
for i=1:size(X,2)
for j=1:size(Y,2)
  Z(j,i)=exp(-X(i)^2-Y(j)^2*0.1);  #行が縦方向Y, 列が横方向X
end
end
#とりあえず描きますよ. で, データを取得するわけですね.
[cnt lev]=contourc(Y,X,Z,levels);
#データを整理します
[val lev]=save_contour(cnt);
#繋いでしまおう
[val lev]=connect_contour(val,lev,1e-5);

で，自作関数 save_contour()とconnect_contour()ですね

function [val lev] = savecontour(cnt)
  npar = 0; ptr = 0; N = size(cnt,2);
  while(ptr < N)
    npar++; ptr++;
    lev{npar} = cnt(1,ptr);
    n = cnt(2,ptr);
    val{npar} = cnt(:,((ptr+1):(ptr+n)))';
    #val{npar}(1:2,1:n)が(X1,Y1)(X2,Y2)...(Xn,Yn)である. だから変換もできる
    #for i=1:n
    #  val{npar}(:,i) =function_x(cnt(:,ptr+i));
    #end みたいな感じで
    ptr+=n;
  endwhile;
endfunction;
function [val lev] = connect_contour(val,lev,errmax)
  len_val=size(val,2);
  j=1;
  do
    k=j+1;
    do
      if (lev{j}==lev{k})
        len_j=size(val{j},2);
        len_k=size(val{k},2);
        if (distance(val{j}(1,len_j)-val{k}(1,len_k),val{j}(2,len_j)-val{k}(2,len_k))<errmax)
          #jの最後とkの最後が繋がるパターン
          for m=len_k:-1:1
            val{j}(:,len_j+len_k-m+1)=val{k}(:,m);
          end
          val(:,k)=[];
          lev(:,k)=[];
          len_val=size(val,2);
          break;
        elseif (distance(val{j}(1,1)-val{k}(1,1),val{j}(2,1)-val{k}(2,1))<errmax)
          #jの最初とkの最初が繋がるパターン
          for m=len_j:-1:1
            val{j}(:,len_k+m)=val{j}(:,m);
          end
          for m=1:len_k
            val{j}(:,m)=val{k}(:,len_k-m+1);
          end
          val(:,k)=[];
          lev(:,k)=[];
          len_val=size(val,2);
          break;
        elseif (distance(val{j}(1,1)-val{k}(1,len_k),val{j}(2,1)-val{k}(2,len_k))<errmax)
          #jの最初とkの最後が繋がるパターン
          for m=len_j:-1:1
            val{j}(:,len_k+m)=val{j}(:,m);
          end
          for m=1:len_k
            val{j}(:,m)=val{k}(:,m);
          end
          val(:,k)=[];
          lev(:,k)=[];
          len_val=size(val,2);
          break;
        elseif (distance(val{j}(1,len_j)-val{k}(1,1),val{j}(2,len_j)-val{k}(2,1))<errmax)
       #jの最後とkの最初が繋がるパターン
         for m=1:len_k
           val{j}(:,len_j+m)=val{k}(:,m);
         end
         val(:,k)=[];
         lev(:,k)=[];
         len_val=size(val,2);
         break;
        endif
      endif
      k++;
    until ( k > len_val);
    j++;
  until ( j > len_val-1);
  #xが増大する向きに揃える
  for j=1:size(val,2)
    len_j=size(val{j},2);
    if (val{j}(1,1)>val{j}(1,len_j))
      for i=1:(len_j/2)
        swp(:,1)=val{j}(:,i);
        val{j}(:,i)=val{j}(:,len_j-i+1);
        val{j}(:,len_j-i+1)=swp(:,1);
      end
    endif
  end
endfunction;