%
% 2次元定在波(TE波)の可視化
% 
% MATLAB R2015a
% 2016/9/14 Copyright (c) Takuichi Hirano
%
function standing_wave_2d_te
    ci=complex(0,1);    % 虚数単位i
    eps0=8.854e-12;     % 真空の誘電率
    mu0=4.*pi.*1e-7;    % 真空の透磁率
    c=2.99792458e8;     % 光速

    lambda0=1.;         % 波長
    freq=c./lambda0;    % 周波数
    omega=2.*pi.*freq;  % 角周波数
    T=1./freq;          % 周期

    %----------------
    % 解析入力パラメータ
    a1=1.;               % 右上に進む波の振幅
    a2=1.;               % 右下に進む波の振幅

    theta=90.*(pi./180.); % 交差角 (rad)
    %----------------
    % 描画・動画パラメータ
    xymin=-1;           % 描画範囲の左・下 (波長単位)
    xymax= 1;           % 描画範囲の右・上 (波長単位)

    n_grid_samp_field=15;    % quiver描画の格子間隔

    canvas_scale=1;     % 描画領域サイズのスケール(1: 画面サイズいっぱい)
    nframe=16;          % 1周期のフレーム数
    nrepeat=10;         % 動画の繰り返し回数

    avi_frate=2;        % フレーム レート (1 秒あたりに表示するフレームの数) 
    avi_quality=70;     % 画質設定 (0 から 100 までのパーセント値)
    avi_nrep=3;         % AVIファイル出力時の動画繰り返し回数
    %----------------

    k0=omega.*sqrt(eps0.*mu0);
    eta0=sqrt(mu0./eps0);

    emax=(abs(a1)+abs(a2));         % 反対方向電界定在波の最大値
    hmax=(abs(a1)+abs(a2))/eta0;    % 反対方向磁界定在波の最大値

    alpha=theta./2;

    % E and H fields
    h1z=@(X,Y) (a1./eta0).*(exp(-ci.*k0.*(cos(alpha).*X+sin(alpha).*Y)));
    h2z=@(X,Y) (a2./eta0).*(exp(-ci.*k0.*(cos(alpha).*X-sin(alpha).*Y)));

    e1x=@(X,Y) a1.*exp(-ci.*k0.*(cos(alpha).*X+sin(alpha).*Y)) ...
        .*(-sin(alpha));
    e1y=@(X,Y) a1.*exp(-ci.*k0.*(cos(alpha).*X+sin(alpha).*Y)) ...
        .*( cos(alpha));
    e2x=@(X,Y) a2.*exp(-ci.*k0.*(cos(alpha).*X-sin(alpha).*Y)) ...
        .*( sin(alpha));
    e2y=@(X,Y) a2.*exp(-ci.*k0.*(cos(alpha).*X-sin(alpha).*Y)) ...
        .*( cos(alpha));

    % グラフィックを描画
    scrsz = get(groot,'ScreenSize');    % 画面サイズの取得
    canvas_width=(scrsz(3))*canvas_scale;
    canvas_height=(scrsz(4))*canvas_scale;
    font_size1=(canvas_width/800)*12;   % フォントサイズ1
    font_size2=(canvas_width/800)*10;   % フォントサイズ2
    fig=figure('Position',[(scrsz(3)/2-canvas_width/2) (scrsz(4)/2-canvas_height/2) canvas_width canvas_height]);

    [X,Y]=meshgrid(xymin:xymax./100:xymax);
    [X2,Y2]=meshgrid(xymin:xymax./n_grid_samp_field:xymax);

    for i=1:nframe
        clf;        % 現在のFigureウィンドウをクリア
        t=(T./nframe).*(i-1);

        % 入射波+反射波
        hz=real((h1z(X,Y)+h2z(X,Y)).*exp(ci.*omega.*t));

        ex=real((e1x(X2,Y2)+e2x(X2,Y2)).*exp(ci.*omega.*t));
        ey=real((e1y(X2,Y2)+e2y(X2,Y2)).*exp(ci.*omega.*t));

        %--------
        % Order of figures
        % 1 2
        %--------
        %---- 1
        subplot(1,2,1)                  % 複数グラフをタイル状に描画
        pcolor(X,Y,hz);                 % グラフを描画
        shading flat;                   % グリッドメッシュを削除
        colormap('jet');                % 色の設定
        caxis([-hmax hmax]);            % カラーマップの範囲を指定
        zlim([-hmax hmax]);             % z軸の範囲を指定
        axis square;                    % グラフ縦横比の設定
        xlabel({'x (wavelength)'});     % xlabel を作成
        ylabel({'y (wavelength)'});     % ylabel を作成
        title({'H-field'});             % title を作成

        hold on;                        % 重ね書きのため
        r=(xymax-xymin)/4;
        v=quiver(0,0,r*cos(alpha),r*sin(alpha),0);           % 平面波方向の矢印を描画
        v.Color='w';
        v.LineWidth=2;
        v.MaxHeadSize=1;
        v=quiver(0,0,r*cos(alpha),-r*sin(alpha),0);           % 平面波方向の矢印を描画
        v.Color='w';
        v.LineWidth=2;
        v.MaxHeadSize=1;
        hold off;
        set(gca,'FontSize',font_size1); % フォントサイズ設定

        %---- 2
        subplot(1,2,2)                  % 複数グラフをタイル状に描画

        hold on;                        % 重ね書きのため
        scale_fact=((3/emax)/(xymax-xymin))/n_grid_samp_field;       % 矢印の長さの自動スケーリング
        h=quiver(X2,Y2,scale_fact*ex,scale_fact*ey);
        h.AutoScale='off';              % 矢印の長さの自動スケーリング
        h.Color='b';                    % 矢印の色
        h.MaxHeadSize=10;
        h.Clipping='on';
        hold off;

        axis square;                    % グラフ縦横比の設定
        axis([xymin xymax xymin xymax]);% データの範囲と描画範囲を一致させる
        xlabel({'x (wavelength)'});     % xlabel を作成
        ylabel({'y (wavelength)'});     % ylabel を作成
        title({'E-field'});             % title を作成

        % 説明用 textbox を描画
        tb1=annotation('textbox',...
            [0.1 0.84 0.35 0.125],...
            'String',{sprintf('Crossing angle (deg): %4.1f',theta.*(180./pi))});
        tb1.FontSize=font_size2;        % フォントサイズ設定
        tb1.EdgeColor='none';           % 枠の色
        %tb1.Color='w';                  % テキストの色
        %tb1.FontWeight='bold';          % 太字
        %tb.BackgroundColor='w';         % 背景色
        %tb.FaceAlpha=0.3;               % 透明度(透明:0～不透明:1)

        set(gca,'FontSize',font_size1); % フォントサイズ設定
        drawnow;                        % 描画する
        frm(i)=getframe(gcf);           % gcf: 現在の Figure ハンドル番号, getframe: ムービーフレームのキャプチャ
    end
    movie(fig,frm,nrepeat);             % 動画の再生(最後の引数は繰り返し回数)

    % 動画AVIファイルの生成
    v = VideoWriter('standing_wave_2d_te.avi');
    v.FrameRate = avi_frate;    % Default 30
    v.Quality = avi_quality;    % Default 75
    open(v);
    for k=1:avi_nrep
        for i=1:nframe
            writeVideo(v,frm(i));
        end
    end
    close(v);
end
