%--------------------------------------------------------------
\documentclass[ja=standard,xelatex,10pt]{bxjsarticle}
%--------------------------------------------------------------

%--------------------------------------------------------------
%packages%{
% #region
\usepackage{bigints}
\usepackage{changepage}
\usepackage[utf8]{inputenc}
\usepackage{lipsum}
\usepackage[bottom]{footmisc}
\usepackage{subcaption}
\usepackage{zxjatype}
\usepackage[hiragino]{zxjafont} 
\usepackage{tcolorbox}
\usepackage{amsmath}
\usepackage{nccmath}
\usepackage{amsthm}
\usepackage{empheq}
\usepackage{amssymb}
\usepackage{bm}
\usepackage{amsfonts}
\usepackage{graphicx}
\usepackage{tocloft}
\usepackage{ascmac}
\usepackage{latexsym}
\usepackage[version=4]{mhchem}
\usepackage{listings, xcolor}
\usepackage{mathtools}
\usepackage[colorlinks=true, linkcolor=blue , citecolor = blue,urlcolor=blue]{hyperref}
\usepackage{array}
\usepackage{url}
\usepackage[english]{babel}
\usepackage[square,numbers]{natbib}
\usepackage{physics}
\usepackage{caption}
\usepackage{ulem}
\usepackage[section]{placeins}
\usepackage{xeCJKfntef}
\usepackage{longtable}
\usepackage{float}
\usepackage{multicol}
\usepackage{mathcomp}
\usepackage{IEEEtrantools}
\usepackage{flafter}
\usepackage{enumitem}
\usepackage[nospread]{cuted}
\setlength{\stripsep}{5pt plus 1pt minus 1pt}
% #endregion
%--------------------------------------------------------------%}

%--------------------------------------------------------------
%define%{
% #region

%%skip,space関係
\newcommand{\vs}[1]{\vskip#1\baselineskip}
\newcommand{\parlskip}[1]{\par\setlength{\leftskip}{#1pt}}

%equation環境
\newcommand{\eqn}[1]{\begin{equation*}{#1}\end{equation*}}
\newcommand{\neqn}[1]{\begin{equation}{#1}\end{equation}}
\newcommand{\leqn}[2]{\begin{fleqn}[#1pt]\begin{equation*}{#2}\end{equation*}\end{fleqn}[#1pt]}
\newcommand{\leqnarray}[2]{\begin{fleqn}[#1pt]\begin{eqnarray*}{#2}\end{eqnarray*}\end{fleqn}}
\newcommand{\eqns}[2]{\begin{equation*}\begin{array}{#1}#2\end{array}\end{equation*}}
\newcommand{\neqns}[2]{\begin{IEEEeqnarray}{#1}#2\end{IEEEeqnarray}}
\newcommand{\neqnarray}[1]{\begin{eqnarray}#1\end{eqnarray}}

%定理・定義・証明等の記述（数学）
\newcommand{\setting}[1]{\vs{0.3}\underline{\textbf{Setting}}\vs{0.3}\parlskip{20}{#1}\parlskip{0}}
\newcommand{\thm}[1]{\vs{0.3}\underline{\textbf{Thm}}\vs{0.3}\parlskip{20}{#1}\parlskip{0}}
\newcommand{\Def}[1]{\vs{0.3}\underline{\textbf{Def}}\vs{0.3}\parlskip{20}{#1}\parlskip{0}}
\newcommand{\rem}[1]{\vs{0.5}\large{\underline{{\textit{remark}}}}\normalsize\vs{0.5}{#1}}
\newcommand{\Proof}[1]{\vs{0.5}\large{\underline{\textit{Proof}}}\normalsize\vs{0.5}{#1}\setcounter{equation}{0}}
\newcommand{\free}[2]{\vs{0.5}\large{\underline{\textit{#1}}}\normalsize\vs{0.5}{#2}\setcounter{equation}{0}}

%tinysections
\newcommand{\tinysection}[3]{\vs{0.5}\large\underline{\textbf{(#1)：#2}}\normalsize\vs{0.5}{#3}}
\newcommand{\tinytinysection}[4]{\vs{0.3}\normalsize{\underline{\textbf{(#1.#2)：#3}}}\vs{0.3}{#4}}
\newcommand{\tinytinytinysection}[5]{\vs{0.2}\normalsize{\underline{(#1.#2.#3)：#4}}\vs{0.2}{#5}}

%微分
\newcommand{\dif}[2]{\dfrac{d#1}{d#2}}
\newcommand{\difop}[1]{\dfrac{d}{d#1}}
\newcommand{\lagdif}[2]{\frac{D#1}{D#2}}

%偏微分
\newcommand{\pdif}[2]{\dfrac{\partial{#1}}{\partial#2}}
\newcommand{\pdifop}[1]{\dfrac{\partial}{\partial#1}}
\newcommand{\pdifl}[1]{\partial_{#1}}
\newcommand{\pdiftwice}[2]{\dfrac{\partial^{2}{#1}}{\partial {#2} ^{2}}}

%積分
\newcommand{\dint}{\displaystyle\int}
\newcommand{\diint}{\displaystyle\iint}
\newcommand{\diiint}{\displaystyle\iiint}

%総和
\newcommand{\sumn}[2]{\sum_{n=#1}^{#2}}
\newcommand{\sumi}[2]{\sum_{i=#1}^{#2}}
\newcommand{\sumj}[2]{\sum_{j=#1}^{#2}}
\newcommand{\sumk}[2]{\sum_{k=#1}^{#2}}
\newcommand{\suml}[2]{\sum_{l=#1}^{#2}}
\newcommand{\summ}[2]{\sum_{m=#1}^{#2}}
\newcommand{\sump}[2]{\sum_{p=#1}^{#2}}
\newcommand{\sumq}[2]{\sum_{q=#1}^{#2}}
\newcommand{\sumfr}[3]{\sum_{#1=#2}^{#3}}
\newcommand{\dsum}{\displaystyle\sum}

%ベクトル
\newcommand{\vect}[1]{\mathbf{#1}}
\newcommand{\vecp}[2]{\mathbf{#1}^{#2}}
\newcommand{\vecn}[1]{\mathbf{#1}^n}
\newcommand{\vecm}[1]{\mathbf{#1}^m}
\newcommand{\pvect}[1]{\begin{pmatrix}#1\end{pmatrix}}

%括弧
\newcommand{\bracketnormal}[1]{\left(#1\right)}
\newcommand{\curlybracket}[1]{\left\{#1\right\}}

%ε-δ論法
\newcommand{\epsilondelta}[4]{\forall\;\epsilon>0,\exists\delta\;\; s.t, \;\; 0<\left| #1-#2 \right|<\delta \Rightarrow \left| #3(#1)-#4\right|<\epsilon}

%その他
\newcommand{\realpower}[1]{\mathbb{R}^{#1}}
\newcommand{\Real}{\mathbb{R}}
\newcommand{\Bf}[1]{\mathbf{#1}}
\newcommand{\limit}[2]{\lim_{#1 \to #2}}
\newcommand{\absdefined}[1]{\left| #1 \right|}
\newcommand{\map}[5]{#1:#2\to#3;\;\;#4\mapsto#5}
\newcommand{\inner}[1]{\left\langle #1 \right\rangle}
\newcommand{\vertline}[2]{\vs{0}\begin{center}\rule{#1\textwidth}{#2pt}\end{center}\vs{0}}
\newcommand{\coord}[2]{\left(#1\;,\;#2\right)}
\newcommand{\scipow}[1]{\times 10^{#1}}
\newcommand{\labeqn}[1]{\label{eqn:#1}}
\newcommand{\labfig}[1]{\label{fig:#1}}
\newcommand{\refeqn}[1]{（\ref{eqn:#1}）}
\newcommand{\reffig}[1]{Figure\ref{fig:#1}}

\usepackage{xcolor}
\usepackage{xparse}
\usepackage{titlesec}

\titlespacing*{\section}{0pt}{1.0\baselineskip}{0.5\baselineskip}
\titlespacing*{\subsection}{0pt}{0.8\baselineskip}{0.4\baselineskip}
\titlespacing*{\subsubsection}{0pt}{0.6\baselineskip}{0.3\baselineskip}



% #endregion
%--------------------------------------------------------------%}

%--------------------------------------------------------------
%title,author,etc%{
% #region
\title{コンピュータグラフィクス論　モデリング課題}
\author{05-242628 三田村彰大（理学部地球惑星物理学科）}
\date{}
% #endregion
%--------------------------------------------------------------%}

\begin{document}

%--------------------------------------------------------------
% #region
%code_area_settings%{
\lstset{
    basicstyle = {\ttfamily\small}, % 基本的なフォントスタイル
    frame = {tbrl}, % 枠線の枠線。t: top, b: bottom, r: right, l: left
    breaklines = true, % 長い行の改行
    numbers = left, % 行番号の表示。left, right, none
    showspaces = false, % スペースの表示
    showstringspaces = false, % 文字列中のスペースの表示
    showtabs = false, %　タブの表示
    keywordstyle = \color{blue}, % キーワードのスタイル。intやwhileなど
    commentstyle = {\color[HTML]{1AB91A}}, % コメントのスタイル
    identifierstyle = \color{black}, % 識別子のスタイル　関数名や変数名
    stringstyle = \color{brown}, % 文字列のスタイル
    captionpos = t % キャプションの位置 t: 上、b: 下
}
\lstdefinelanguage{JavaScript}{
  keywords={const, let, var, function, return, if, else, for, while, do, switch, case, break, continue, new, class, import, from, export, default, null, true, false},
  keywordstyle=\color{blue},
  ndkeywords={document, window, console, Math, Number, Array, Object, THREE},
  ndkeywordstyle=\color{teal},
  identifierstyle=\color{black},
  sensitive=true,
  comment=[l]{//},
  morecomment=[s]{/*}{*/},
  commentstyle=\color[HTML]{1AB91A},
  stringstyle=\color{brown},
  morestring=[b]',
  morestring=[b]",
  morestring=[b]`
}
\captionsetup*[lstlisting]{labelformat=empty}
% #endregion
%--------------------------------------------------------------%}

%--------------------------------------------------------------
%document_settings%{
% #region

\captionsetup*[figure]{labelfont=bf,font=footnotesize}
\captionsetup*[table]{labelfont=bf,font=footnotesize}
\captionsetup*[subfigure]{position=bottom,font=footnotesize}
\renewcommand{\arraystretch}{1.15}

\setlength{\dbltextfloatsep}{10pt plus 2pt minus 2pt}
\setlength{\dblfloatsep}{8pt plus 2pt minus 2pt}

\setcounter{tocdepth}{3}
\setlength{\parindent}{0pt}
\newcommand{\setdisplaymathskip}{%
\setlength{\abovedisplayskip}{10pt plus 1pt minus 2pt}
\setlength{\belowdisplayskip}{10pt plus 1pt minus 2pt}
\setlength{\abovedisplayshortskip}{6pt plus 1pt minus 1pt}
\setlength{\belowdisplayshortskip}{6pt plus 1pt minus 1pt}
\setlength{\jot}{2pt}
}
\setlength{\jot}{3pt}
\newgeometry{top=18truemm,bottom=20truemm,left=14truemm,right=14truemm}
\renewcommand{\arraystretch}{1.5}
\bibliographystyle{abbrvnat}
\setlength{\columnsep}{5mm}
\setlength{\columnseprule}{0pt}
%asm_settings
\theoremstyle{definition}
\newtheorem{prop}{proposition}[section]
%subfigure_settings
\captionsetup*[subfigure]{position=bottom}
\captionsetup*[figure]{labelfont=bf}
\captionsetup*[table]{labelfont=bf}
%footnote_settings
\setlength{\footnotesep}{0.4\baselineskip} % 脚注どうしの間隔
\makeatletter
\renewcommand\@makefntext[1]{%
  % ---- 脚注本文の段落設定（同一脚注内） ----
  \setlength{\parindent}{0pt}%
  \setlength{\parskip}{0.4\baselineskip}%
  % ---- ぶら下げ＋マーク表示 ----
  \noindent
  \hb@xt@1.8em{\hss\@makefnmark}%
  \hspace{0.4em}%
  \hangindent=2.2em%
  \hangafter=1%
  #1%
}
\renewcommand{\footnoterule}{%
  \kern -3pt
  \hrule width 1.0\linewidth
  \kern 5pt
}
\makeatother
\makeatletter
\renewcommand{\bibsection}{} %"References" 見出しを無効化
\makeatother
\makeatletter
\renewenvironment{proof}[1][\proofname]{%
  \par\pushQED{\qed}\normalfont
  \trivlist
  \item[\hskip\labelsep\bfseries #1\@addpunct{.}]%
  \ignorespaces
}{%
  \popQED\endtrivlist\@endpefalse
}
\makeatother
% #endregion
\renewcommand{\abstractname}{}


%--------------------------------------------------------------%}


%--------------------------------------------------------------
%maketitle_makehedding_tableofcontents%{
% #region

\twocolumn[{\maketitle


\begin{description}[
  labelwidth=4.5cm,
  labelsep=0.5cm,
  leftmargin=5.0cm,
  style=multiline
]
  \item[動作が確認できるURL：] \href{https://shota-mitamura.com/works/computergraphicskadai1/template}{https://shota-mitamura.com/works/computergraphicskadai1/template}
  \item[コード置き場：] \href{https://github.com/ramutami/modeling_mitamura}{https://github.com/ramutami/modeling\_mitamura}
  \item[AIとのチャットログ：] \href{https://chatgpt.com/g/g-p-69fabe32fd9c81919d69b22d4294ca7e-konhiyutakurahuikusulun-ke-ti/project}{https://chatgpt.com/g/g-p-69fa...}
\end{description}

\vs{2}
}]



\setdisplaymathskip



% #endregion
%--------------------------------------------------------------%}

%start_document



\section{概要}
　パラメトリック曲線、特にCatmull-Romスプラインを、パラメタがUniformな場合とCentripetalな場合で実装し、両者を比較する。特に、Catmull-Romスプラインが制御点を滑らかに補間することや、Centripetalなパラメタ区分を用いることで望ましくない挙動を抑制できることなどを確認する。
\subsection*{※AIの利用について}
　本レポートでは、パラメトリック曲線の実装にあたり.html/.css/.jsファイルの雛形をAIに相談しながら作成した\footnote{
  普段FORTRANしか触っていないので…
}。具体的には、動作確認用のWebsiteの見た目や、制御点のドラッグ操作、スライダによる制御点の個数変更などのUIの実装周りをChatGPTを用いたバイブコーディングで行った。その上で実際のパラメトリック曲線の実装やパラメタ区分の実装は自分で行うこととした\footnote{
  AIの察しが良すぎたがために、UIの相談をしているだけであるにもかかわらずパラメトリック曲線の実装コードをサジェストされてしまったが、できるだけそれを無視してコーディングを行った。
}。


\section{理論}
\subsection{Catmull-Romスプライン}
　時刻$t_{k}$での値が$\vect{x}_k$で表されるような状況を考える。（\reffig{1}参照）この時、$\vect{x}_{k}$と$\vect{x}_{k+1}$を繋ぐ曲線$r_k(t)$（図中では$x_k(t)$として表示）を$\vect{x}_{k-1},\vect{x}_{k},\vect{x}_{k+1},\vect{x}_{k+2}$を用いて次のように決定する。
\begin{enumerate}
  \item 直線による補間を行う。（\reffig{2}参照）
  \begin{itemize}
    \item $\vect{x}_{k-1}\sim \vect{x}_k$を結ぶ直線$l_{k-1}(t)$を求める。
    \item  $\vect{x}_{k}\sim \vect{x}_{k+1}$を結ぶ直線$l_{k}(t)$を求める。
    \item $\vect{x}_{k+1}\sim \vect{x}_{k+2}$を結ぶ直線$l_{k+1}(t)$を求める。
  \end{itemize}
  \item 二次曲線による補間を行う。（\reffig{3}参照）
  \begin{itemize}
    \item $l_{k-1}$と$l_{k}$を用いて、$\vect{x}_{k-1}\sim \vect{x}_{k}\sim \vect{x}_{k+1}$の間を二次曲線$q_k(t)$で補間する。
    \item $l_{k}$と$l_{k+1}$を用いて、$\vect{x}_{k}\sim \vect{x}_{k+1}\sim \vect{x}_{k+2}$の間を二次曲線$q_{k+1}(t)$で補間する。
  \end{itemize}
  \item 三次曲線による補間を行う。（\reffig{4}参照）
  \begin{itemize}
    \item $q_k$と$q_{k+1}$を用いて、$\vect{x}_k\sim \vect{x}_{k+1}$の間を三次曲線$r_k(t)$で補間する。
  \end{itemize}
\end{enumerate}

　ただし、それぞれの曲線（直線）の式は次のように定める。
\neqns{lcl}{
  l_k(t) &\; =\; & \bracketnormal{1-\dfrac{t-t_k}{t_{k+1}-t_{k}}}\vect{x}_k + \dfrac{t-t_k}{t_{k+1}-t_k}\vect{x}_{k+1} \labeqn{cat1}\\[10pt]
  q_k(t) &\; =\; & \bracketnormal{1-\dfrac{t-t_{k-1}}{t_{k+1}-t_{k-1}}}l_{k-1}(t) \nonumber\\
  && \qquad \qquad \qquad \qquad  + \dfrac{t-t_{k-1}}{t_{k+1}-t_{k-1}}l_k(t)\labeqn{cat2} \\[10pt]
  r_k(t) &\; =\; &\bracketnormal{1-\dfrac{t-t_k}{t_{k+1}-t_k}}q_k(t) \nonumber \\
  && \qquad \qquad \qquad \qquad+ \dfrac{t-t_k}{t_{k+1}-t_k}q_{k+1}(t)\labeqn{cat3}
}
　よって実際に実装する上では、$\vect{x}_{k-1},\vect{x}_{k},\vect{x}_{k+1},\vect{x}_{k+2}$及び$t_{k-1},t_{k},t_{k+1},t_{k+2}$及び$t$を入力するとパラメトリック曲線上の座標$r_k(t)$を返すような次のような関数を定義する\footnote{
  描画は二次元空間上に行うが、実装上は$z=0$とした三次元空間内で描画するので、実際の関数としてはThree.Vector3(x,y,z)を返すものとなる。（一応コードでは、今後の拡張も考えて$z$座標の更新もできるようになっているが、結局使わなかった…）
}。

\small
\neqns{l}{
  \text{function}(t , \vect{x}_{k-1},\vect{x}_{k},\vect{x}_{k+1},\vect{x}_{k+2},t_{k-1},t_{k},t_{k+1},t_{k+2})\{\nonumber \\
  \quad \text{calculate:}\; l_{k_i}\;\text{using}\; t,t_{k_i},\vect{x}_{k_i}\nonumber \\
  \quad \text{calculate:}\; q_{k_i}\;\text{using}\; t,t_{k_i},l_{k_i}\nonumber \\
  \quad \text{calculate:}\; r_{k}\;\text{using}\; t,t_{k_i},q_{k_i}\nonumber \\
  \quad \text{return:}\; r_k(t \, ; \, \vect{x}_{k-1},\vect{x}_{k},\vect{x}_{k+1},\vect{x}_{k+2},t_{k-1},t_{k},t_{k+1},t_{k+2}) \nonumber\\
  \text{\}}\labeqn{func}
}
\normalsize
\setdisplaymathskip


\begin{figure}[!htbp]
  \centering
  \includegraphics[width=\columnwidth]{catmull.png}
  \caption{考えている問題設定。授業スライドより抜粋。}
  \label{fig:1}

  \includegraphics[width=0.8\columnwidth]{ichiji.png}
  \caption{直線による補間。授業スライドより抜粋。}
  \label{fig:2}

  \includegraphics[width=0.8\columnwidth]{niji.png}
  \caption{二次曲線による補間。授業スライドより抜粋。}
  \label{fig:3}

  \includegraphics[width=0.8\columnwidth]{sanji.png}
  \caption{三次曲線による補間。ただし文中では$x_k(t)$ではなく$r_k(t)$としている。授業スライドより抜粋。}
  \label{fig:4}
\end{figure}

\FloatBarrier
\section{実装}
\subsection{Uniformなパラメタ区分}
　$t_k = t_{k-1}+1$で表されるようなUniformなパラメタ区分を考える。このとき、一般に$t_k=0$としても問題なく、
\neqn{
  t_{k-1}=-1\; ,\;t_k=0\; ,\;t_{k+1}=1\; ,\;t_{k+2}=2
}
を\refeqn{cat1}〜\refeqn{cat3}に代入しながら先のアルゴリズムに従って$r_k(t)$を決定すると最終的なパラメタ曲線の座標は次のように具体的に計算することができる。
\vs{0}
\footnotesize
\neqns{l}{
  r_k(t \, ; \, \vect{x}_{k-1},\vect{x}_{k},\vect{x}_{k+1},\vect{x}_{k+2})\nonumber \\[10pt]
  \quad = \bracketnormal{1-t}\bigg[\bracketnormal{1-\dfrac{t+1}{2}}\bigg\{\bracketnormal{1-\bracketnormal{t+1}}\vect{x}_{k-1} + \bracketnormal{t+1}\vect{x}_k\bigg\}\nonumber \\
  \quad\qquad\qquad\qquad + \bracketnormal{\dfrac{t+1}{2}}\bigg\{\bracketnormal{1-t}\vect{x}_k + t \;\vect{x}_{k+1}\bigg\}\bigg]\nonumber\\
  \quad \quad + t \;\bigg[\bracketnormal{1-\dfrac{t}{2}}\bigg\{ \bracketnormal{1-t}\vect{x}_k + t \;\vect{x}_{k+1} \bigg\} \nonumber \\
  \quad\qquad\qquad\qquad + \dfrac{t}{2}\;\bigg\{ \bracketnormal{1-\bracketnormal{t-1}}\vect{x}_{k+1} + \bracketnormal{t-1}\vect{x}_{k+2} \}\bigg]
}
\normalsize
\setdisplaymathskip
\vs{0}
　これはさらに整理することができて、\href{https://ja.wikipedia.org/wiki/Catmull-Rom\%E3\%82\%B9\%E3\%83\%97\%E3\%83\%A9\%E3\%82\%A4\%E3\%83\%B3\%E6\%9B\%B2\%E7\%B7\%9A}{Wikipedia}に記載されているものと同様の次の式を得ることができる。

\small
\neqn{
  r_k(t) = \dfrac{1}{2}\begin{bmatrix}
    t^3 & t^2 & t & 1
  \end{bmatrix}\begin{bmatrix}
    -1 & 3 & -3 & 1 \\
    2 & -5 & 4 & -1 \\
    -1 & 0 & 1 & 0 \\
    0 & 2 & 0 & 0
  \end{bmatrix}\begin{bmatrix}
    \vect{x}_{k-1}\\
    \vect{x}_{k}\\
    \vect{x}_{k+1}\\
    \vect{x}_{k+2}\\
  \end{bmatrix}
}
\normalsize
\setdisplaymathskip
\vs{0}
　これを実装したのが、\href{https://github.com/ramutami/modeling_mitamura/blob/main/catmuluniform.js}{catmuluniform.js} 内の関数 catmullRomuniform である。\footnote{
  コード57行目。
}
\begin{lstlisting}[
  language=JavaScript,
  numbers=none,
  title={catmullRomuniform関数},
  caption={catmullRomuniform関数},
  basicstyle=\ttfamily\footnotesize,
  xleftmargin=1em
]
function catmullRomuniform(pkm1, pk, pkp1, pkp2, t) {
  const t2 = t * t;
  const t3 = t2 * t;

  const x = 0.5*(
    t3*(-1*pkm1.x+3*pk.x-3*pkp1.x+pkp2.x)+
    t2*(2*pkm1.x-5*pk.x+4*pkp1.x-pkp2.x)+
    t*(-pkm1.x+pkp1.x)+
    2*pk.x
  );

  const y = ...
  const z = ...

  return new THREE.Vector3(x, y, z);
}
\end{lstlisting}

　よって points[ ] で表されるような制御点$\vect{x}_0,\vect{x}_1,\ldots$に対し、実際にパラメトリック曲線を描画する部分は次のようすれば良いとわかる。

\begin{lstlisting}[
  language=JavaScript,
  numbers=none,
  title={パラメトリック曲線の描画},
  caption={パラメトリック曲線の描画部分},
  basicstyle=\ttfamily\footnotesize,
  xleftmargin=1em
]
for (let i = 0; i < points.length - 1; i++) {
  const pk_1 = points[i-1];
  const pk = points[i];
  const pkp1 = points[i+1];
  const pkp2 = points[i+2];
  for (let j = 0; j < samplesPerSegment; j++) {
    const t = j / samplesPerSegment;
    curvePoints.push(catmullRomuniform(pk_1, pk, pkp1, pkp2, t));
  }
}
\end{lstlisting}

　すなわち制御点$\vect{x}_{i-1}\sim \vect{x}_{i+2}$を用いてパラメトリック曲線を計算し、時刻$t=0+j\Delta t$での曲線の座標を計算すれば良い。\footnote{
  $t_i = 0$としたことを利用している。
}

\subsection{端点処理}
　最初の制御点$\vect{x}_0$と$\vect{x}_1$の間のパラメトリック曲線を描画する際、すなわち$\vect{x}_k = \vect{x}_0$の時は、$\vect{x}_{k-1}$として用いることのできる制御点が存在しない。そこでこの場合は、$\vect{x}_{k-1} = \vect{x}_{k}=\vect{x}_0$として処理することにする。\footnote{
  最後の制御点$\vect{x}_{N-1}$と$\vect{x}_N$の間のパラメトリック曲線の描画についても同様に処理する。
}
\vs{0}
　これはすなわち、$\vect{x}_0$と同じ値を持つ仮想的な制御点$\vect{x}_{-1}$を追加していることに他ならず、以降の計算アルゴリズムは変更を受けない。（\reffig{tantenuni}参照。）（強いて言えば$\bracketnormal{t,\vect{x}}$空間で$l_{k-1}$直線の傾きが平坦になるだけである。）

\begin{figure}[!htbp]
  \centering
  \includegraphics[width=\columnwidth]{uni_tanten.png}
  \caption{ユニフォームなパラメタ区分における端点処理。}
  \label{fig:tantenuni}
\end{figure}
\FloatBarrier

\subsection{実装結果（Uniformなパラメタ区分）}

　ユニフォームなパラメタ区分を用いた場合のCatmull-Rom曲線を描画した図が\reffig{catuni1}及び\reffig{catuni2}である。まず\reffig{catuni1}を見ると、全体的に滑らか（$C^1$連続）な曲線となっており、端点付近でもその滑らかさが失われていないことがわかる。また、図中央の制御点が非常に密になっている箇所では、「くるりん」とでも言うべき望ましくない挙動（ループ）が現れている。この挙動は後ほどCentripetalなパラメタ区分で解消されることが期待される。さらに、図右側中央の制御点が四角形を形成している箇所に着目すると、この場所でパラメトリック曲線は「ぼわっと」膨らんだ四角形を描いている。このような挙動を見ると、Catmull-Romスプラインがシャープな角の表現を苦手としていることが推察される。
\vs{0}
　また、\reffig{catuni2}はネコの絵をCatmull-Rom曲線で描いたものである。制御点をパラメトリック曲線が必ず通るため、非常に直感的かつ簡単な操作で絵を描くことができた。ただしこの場合においても、$C^1$連続性を保つと言う要請が逆にシャープな角の描画を困難にしており、ネコの「耳」の接続部で不自然な凹みが生じてしまっていることがわかる。この点はcentripetalで改善されるかどうかも注目する必要がある。


\begin{figure}[!htbp]
  \centering
  \includegraphics[width=0.8\columnwidth]{uniform_pic1.png}
  \caption{ユニフォームなパラメタ区分を用いたCatmull-Rom曲線。}
  \label{fig:catuni1}
  \includegraphics[width=0.8\columnwidth]{neko_uniform.png}
  \caption{ネコ（？）の絵。}
  \label{fig:catuni2}
\end{figure}

\FloatBarrier

\subsection{Centripetalなパラメタ区分}
　次に、centripetalなパラメタ区分、すなわち$t_k = t_{k-1} + \sqrt{\norm{\vect{x}_{k-1}-\vect{x}_k}}$と表されるような場合にCatmull-Rom曲線を実装することを考える。この場合、Uniformなパラメタ区分で行ったように$t_{k-1},t_k,t_{k+1},t_{k+2}=-1,0,1,2$を用いて計算を簡略化することはできず、愚直に\refeqn{func}で示したような逐次的に$r_k(t)$を計算する関数を定義する必要がある。（ただし、$t_k=0$は用いることにする。）それを実装したのが \href{https://github.com/ramutami/modeling_mitamura/blob/main/catmulcentripetal.js}{catmulcentripetal.js} 内の catmullRomcentripetal 関数\footnote{コード58行目}である。


\begin{lstlisting}[
  language=JavaScript,
  numbers=none,
  title={catmullRomcentripetal},
  caption={catmullRomcentripetal関数},
  basicstyle=\ttfamily\footnotesize,
  xleftmargin=1em
]
function catmullRomcentripetal(pkm1, pk, pkp1, pkp2, tkm1,tkp1,tkp2,t) {
    // lkの計算(x座標)
    const lkm1x = (1-(t-tkm1)/(-tkm1))*pkm1.x + ((t-tkm1)/(-tkm1))*pk.x;
    const lkx = (1-t/tkp1)*pk.x + (t/tkp1)*pkp1.x;
    const lkp1x = (1-(t-tkp1)/(tkp2-tkp1))*pkp1.x + ((t-tkp1)/(tkp2-tkp1))*pkp2.x;
    //qkの計算(x座標)
    const qkx = (1-(t-tkm1)/(tkp1-tkm1))*lkm1x+ ((t-tkm1)/(tkp1-tkm1))*lkx;
    const qkp1x = (1-t/tkp2)*lkx+ (t/tkp2)*lkp1x;
    //rkの計算(x座標)
    const x = (1-t/tkp1)*qkx + (t/tkp1)*qkp1x;
    ... 
    const y = ...
    const z = ...

  return new THREE.Vector3(x, y, z);
}
\end{lstlisting}

　次に、この catmullRomcentripetal 関数を用いて実際にパラメトリック曲線を描画していくわけであるが、その際制御点の座標から$t_{k-1},t_{k+1},t_{k+2}$の値を次のように計算する必要がある。\footnote{
  $t_k = 0$としていることに注意。
}
\neqns{l}{
  t_{k-1} = 0 - \sqrt{\norm{\vect{x}_{k-1} - \vect{x}_{k}}}\\[3pt]
  t_{k+1} = 0 + \sqrt{\norm{\vect{x}_k - \vect{x}_{k+1}}}\\[3pt]
  t_{k+2} = t_{k+1} + \sqrt{\norm{\vect{x}_{k+1} - \vect{x}_{k+2}}}
}
　また、実際にパラメトリック曲線を描画する上では、時刻$t$のグリッド幅$\Delta t$をあらかじめ定めておき、$t_k=0\sim t_{k+1}$の間を$t_{k+1}/\Delta t$個のメッシュに区切って描画する。こうすることで、制御点間の距離が短い時は柔軟にパラメタの間隔を小さくしループなどの望ましくない挙動を抑えることができると期待される。
\vs{0}
　すなわち、実際にパラメトリック曲線を描画する部分は次のようになる。

\small
\neqns{l}{
  \text{loop for }k\text{:}\;\{\nonumber \\
  \quad \text{calculate:}\; t_{k-1} = 0 - \sqrt{\norm{\vect{x}_{k-1} - \vect{x}_{k}}}\nonumber \\
  \quad \text{calculate:}\; t_{k+1} = 0 + \sqrt{\norm{\vect{x}_k - \vect{x}_{k+1}}}\nonumber \\
  \quad \text{calculate:}\; t_{k+2} = t_{k+1} + \sqrt{\norm{\vect{x}_{k+1} - \vect{x}_{k+2}}}\nonumber \\
  \quad \text{loop for }t \in \left[ 0,\Delta t, 2\Delta t,\ldots ,t_{k+1}\right]\{\nonumber\\
  \quad\quad \text{add new point on line} \nonumber\\
  \quad\quad \quad \leftarrow r_k(t \, ; \, \vect{x}_{k-1},\vect{x}_{k},\vect{x}_{k+1},\vect{x}_{k+2},t_{k-1},t_k=0,t_{k+1},t_{k+2}) \nonumber\\
  \quad \}\nonumber\\
  \}\labeqn{func2}
}
\normalsize
（ただし$r_k$として先に示した catmullRomcentripetal 関数を用いることになる。）
\vs{0}
　これを実際に実装すると次のよう。

\begin{lstlisting}[
  language=JavaScript,
  numbers=none,
  title={catmullRomcentripetal},
  caption={パラメトリック曲線(centripetal)の描画部分},
  basicstyle=\ttfamily\footnotesize,
  xleftmargin=1em
]
for (let i = 0; i < points.length - 1; i++) {
  const pkm1 = points[i-1];
  const pk = points[i];
  const pkp1 = points[i+1];
  const pkp2 = points[i+2];

  const tkm1 = -Math.sqrt(Math.sqrt((pkm1.x-pk.x)**2 + (pkm1.y-pk.y)**2 + (pkm1.z-pk.z)**2 ))
  const tkp1 = Math.sqrt(Math.sqrt((pkp1.x-pk.x)**2 + (pkp1.y-pk.y)**2 + (pkp1.z-pk.z)**2 ))
  const tkp2 = tkp1 + Math.sqrt(Math.sqrt((pkp1.x-pkp2.x)**2 + (pkp1.y-pkp2.y)**2 + (pkp1.z-pkp2.z)**2 ))

  const samplesPerSegment = tkp1/deltat
  for (let j = 0; j < samplesPerSegment; j++) {
    const t = j*deltat;
    curvePoints.push(catmullRomcentripetal(pkm1, pk, pkp1, pkp2, tkm1,tkp1,tkp2,t));
  }
}
\end{lstlisting}

\subsection{端点処理}
　Centripetalなパラメタ区分の場合においても、基本的にはUniformなパラメタ区分の場合と同様の端点処理をすることになる。すなわち、最初の制御点$\vect{x}_0$と$\vect{x}_1$の間のパラメトリック曲線を描画する場合は$\vect{x}_{k-1}=\vect{x}_k=\vect{x}_0$として描画するわけであるが、Centripetalな場合においてはこの時$t_{k-1}=t_k=0$となってしまうことに注意しなければならない。というのも、$r_k(t)$を逐次的に計算する際、$l_{k-1}$の計算で$t_k-t_{k-1}$で割ると言う操作が出てきてしまうため$t_{k-1}=t_k$となると都合が悪い。（下式参照）
\neqn{
  l_{k-1}(t)  = \bracketnormal{1-\dfrac{t-t_{k-1}}{t_{k}-t_{k-1}}}\vect{x}_{k-1} + \dfrac{t-t_{k-1}}{t_{k}-t_{k-1}}\vect{x}_{k}
}
　そこで、$k=0$の場合は$l_{k-1}(t) = \vect{x}_{k-1}$として直線$l_{k-1}$を定数関数で表すことにする。（\reffig{catuni25}参照）こうすれば$l_{k-1}$と$l_k$の補間$q_k$はこれまで通り二次曲線となり、以降のアルゴリズムは特別変更を受けない。（\reffig{catuni3}参照）\footnote{
  最後の制御点$\vect{x}_{N-1}$と$\vect{x}_N$の間のパラメトリック曲線の描画についても同様に処理する。
}

\begin{figure}[!htbp]
  \centering
  \includegraphics[width=0.85\columnwidth]{tanten3.pnb.png}
  \caption{Centripetalなパラメタ区分における端点処理。}
  \label{fig:catuni25}

  \includegraphics[width=0.75\columnwidth]{tanten2.png}
  \caption{端点処理後の二次関数による補間。}
  \label{fig:catuni3}
\end{figure}

\subsection{実装結果（Centripetalなパラメタ区分）}

\begin{figure}[!htbp]
  \centering
  \includegraphics[width=0.8\columnwidth]{centri.png}
  \caption{Centripetalなパラメタ区分を用いたCatmull-Rom曲線。}
  \label{fig:centri1}
  \includegraphics[width=0.7\columnwidth]{centrineko.png}
  \caption{Centripetalなパラメタ区分を用いた場合のネコ（？）の絵。}
  \label{fig:catcentri2}
\end{figure}
　Centripetalなパラメタ区分を用いた場合のCatmull-Rom曲線を描画した図が\reffig{centri1}と\reffig{catcentri2}である。まず\reffig{centri1}をみると、大まかな挙動はパラメタがユニフォームな場合（\reffig{catuni1}）と変わっていないが、例えば図中央の制御点が密になっている箇所でループが生じていないなどの変化があることがわかる。すなわち、Centripetalなパラメタ区分を用いたことで望ましくない挙動を抑えることができている。
\vs{0}
　ただし、\reffig{centri1}右側の制御点が四角形をなしている箇所や、左側のギザギザした場所を見ると、ユニフォームなパラメタ区分を用いた場合と同様「ぼわっと」曲線が膨らんでしまっていることがわかる。この点についてはのちの章で考察する。
\vs{0}
　また\reffig{catcentri2}はネコの絵をCentripetalなパラメタ区分のCatmull-Rom曲線で描いた図である。\reffig{catuni2}と比べると、ネコの「耳」の接続部に生じていた不自然な凹みがそれなりに解消されていることがわかる。ただし実際に絵を描いた感想としては、ユニフォームなパラメタ区分を用いた場合に比べて曲線が「ぼわっと」広がってしまうために望ましい曲線を得るのが困難だと感じた。
\section{手法の比較}
\subsection{遠心力的な膨らみの抑制}
　Catmull-Rom曲線においては、曲線の折り曲がりに際してその折れ曲がりとは逆方向へ曲線が膨らむ「遠心力」のような現象が見られることがわかった。
\vs{0}
　これを示したのが\reffig{blow}である。今ユニフォームなパラメタ区分で曲線を描いた\reffig{blow}左図を見ると、左から三番目の制御点の手前で曲線が折れ曲がりの向きとは逆方向に膨らんでしまっていることがわかる。このような現象はシャープな表現を困難にすると思われ、解消されることが望ましい。実際、Centripetalなパラメタ区分を用いた場合である\reffig{blow}右図を見ると、「膨らみ」が抑えられていることがわかりCentripetalなパラメタ区分の有用性が見て取れる。

\begin{figure}[!htbp]
    \centering
    \begin{subfigure}{0.4\columnwidth}
        \centering
        \includegraphics[width = \linewidth]{uni4.png}
    \end{subfigure}
    \hspace{0.05\columnwidth}
    \begin{subfigure}{0.45\columnwidth}
        \centering
        \includegraphics[width=\linewidth]{centri4.png}
    \end{subfigure}
    \caption{折れ曲がりに際して生じる折れ曲がりとは逆方向への曲線の膨らみ。左図：Uniformなパラメタ区分。右図：Centripetalなパラメタ区分}
    \label{fig:blow}
\end{figure}
\FloatBarrier

\subsection{望ましくないループの抑制}
　次に、Centripetalなパラメタ区分によるループの抑制について考察する。ユニフォームなパラメタ区分において、制御点が密になっている箇所で望ましくないループが生じることはこれまで見てきた通りであり、それを示したのが\reffig{looploop}左図である。（図中央は制御点が二つほぼ重なっている。）これは制御点間の距離が短いにも関わらず時間間隔$t$が短くなっていないために無理やり曲線を描こうとしてしまうために生じる現象であり、Centripetalなパラメタ区分を用いることで抑制することができる。（\reffig{looploop}右図参照。見にくいが、図中央は二つの制御点がほとんど重なっているにも関わらずループが生じていない。）

\begin{figure}[!htbp]
    \centering
    \begin{subfigure}{0.4\columnwidth}
        \centering
        \includegraphics[width = \linewidth]{uni5.png}
    \end{subfigure}
    \hspace{0.05\columnwidth}
    \begin{subfigure}{0.45\columnwidth}
        \centering
        \includegraphics[width=\linewidth]{centri5.png}
    \end{subfigure}
    \caption{制御点が密に重なっている時のCatmull-Rom曲線の挙動。左図：Uniformなパラメタ区分。右図：Centripetalなパラメタ区分}
    \label{fig:looploop}
\end{figure}
\FloatBarrier

\subsection{曲がり角が連続する領域}
　\reffig{square}は制御点が四角形を成している時のCatmull-Rom曲線を描いたものである。このような$C^1$不連続な曲がり角が連続する場合は、Catmull-Rom曲線の$C^1$連続性が逆に仇となり曲線が「ぼわっと」膨らんでしまうことがわかる。また、これはパラメタの取り方の問題ではなく$C^1$連続性というCatmull-Rom曲線の根本設計によって生じる問題であるため、Centripetalなパラメタ区分を採用しても解消されない。（\reffig{square}右図参照。）
\vs{0}
　そこで今回は、制御点を密にとることによりこの問題を解消することを試みた。これを表したのが\reffig{sqsq}である。この図においては、曲がり角付近で制御点を密にすることでシャープな折れ曲がりの連続を表現できるようになっている。
\vs{0}
　ただしUniformなパラメタ区分においては、制御点間の距離が大きくなると曲がり角付近の制御点を密にしたことによってループの問題が出てきてしまう点に注意する必要がある。（\reffig{unisqsqsq}参照。）


\begin{figure}[!htbp]
    \centering
    \begin{subfigure}{0.3\columnwidth}
        \centering
        \includegraphics[width = \linewidth]{uni_square.png}
    \end{subfigure}
    \hspace{0.05\columnwidth}
    \begin{subfigure}{0.35\columnwidth}
        \centering
        \includegraphics[width=\linewidth]{centrisquare.png}
    \end{subfigure}
    \caption{制御点が四角形を成している時のCatmull-Rom曲線の挙動。左図：Uniformなパラメタ区分。右図：Centripetalなパラメタ区分}
    \label{fig:square}
\end{figure}

\begin{figure}[!htbp]
    \centering
    \begin{subfigure}{0.35\columnwidth}
        \centering
        \includegraphics[width = \linewidth]{unisqsq.png}
    \end{subfigure}
    \hspace{0.05\columnwidth}
    \begin{subfigure}{0.3\columnwidth}
        \centering
        \includegraphics[width=\linewidth]{centrisqsq.png}
    \end{subfigure}
    \caption{制御点を密にすることでシャープな曲がり角の連続を描くことができる様子。左図：Uniformなパラメタ区分。右図：Centripetalなパラメタ区分}
    \label{fig:sqsq}
\end{figure}


\begin{figure}[!htbp]
  \centering
  \includegraphics[width=0.4\columnwidth]{unisqsqsq.png}
  \caption{制御点間隔が大きい場合に、シャープな折れ曲がりの連続を折れ曲がり付近の制御点を密にすることで表現しようとした図。パラメタ区分はUniform。この図を見ると、先に述べたループ現象の影響により折れ曲がりの点でGibbsの角のようなものが現れてしまっている。また、右下の角は曲がり角付近の制御点間隔を少し広く取ることでループを抑えようとしているが、今度は逆にシャープな角の表現が難しくなってしまう。}
  \label{fig:unisqsqsq}


\end{figure}
\FloatBarrier

\section{おわりに}
　今回はCatmull-Rom曲線を実装し、その挙動を確認した。その中で、Catmull-Rom曲線が実際に滑らかな曲線となっていることや、制御点を必ず通るという制約により非常に直感的かつ簡単に絵を描くことができることなどがわかった。
\vs{0}
　また、Uniformなパラメタ区分とCentripetalなパラメタ区分の比較では、後者を用いることにより前者で現れる望ましくない挙動を抑制できることがわかった。ただし、Catmull-Rom曲線の$C^1$連続性に起因するシャープな曲がり角の描画の困難さはCentripetalなパラメタ区分でも解消できないこともわかった。

\section*{参考文献}
[1] ja.Wikipedia Catmull-Romスプライン曲線 \url{https://ja.wikipedia.org/wiki/Catmull-Rom%E3%82%B9%E3%83%97%E3%83%A9%E3%82%A4%E3%83%B3%E6%9B%B2%E7%B7%9A}




\end{document}



