Informationsverarbeitung in Sensornetzwerken
Prof. Dr.-Ing. Uwe D. Hanebeck
%
Dipl.-Ing. Frederik Beutler
Dipl.-Ing. Florian Wei�el
In diesem Dokument werden die Macros erklärt, welche in den Dateien ISASMatheMacros.sty und DietrichsMatheMacros.sty definiert sind. Sie enthalten Befehle, welche der am ISAS gebräuchlichen Konventionen für Formeln und Formelzeichen entsprechen.
ISASPubl_Keys.tex
\renewcommand{\refname}{\normalsize \Deutsch{Preprints}\English{Preprints} (6)}
\begin{bibunit}
  \nocite{
  	% From here 06/24/12
	arXiv12_Hekler,
	arXiv12_Hanebeck,
	arXiv12_Deisenroth,
	arXiv11_Kretz,
	arXiv11_Lyons,
	arXiv10_Perez}
  	% To here 06/24/12
  \putbib
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Artikel in Konferenzb\"anden}\English{Conference Papers} (218)}
\begin{bibunit}
  \nocite{
  	% From here 11/02/12
  	CDC12_Baum,
	CDC12_Chlebek,
	CDC12_Hekler,
	CDC12_Noack,
	CDC12_Reinhardt,
	IROS12_Faion,
	MFI12_Baum,
	MFI12_Reinhardt,
  	% To here 11/02/12  
  	% From here 06/24/12
	Fusion12_Baum,
	Fusion12_BaumWillett,
	Fusion12_Faion-CylinderTracking,
	Fusion12_Faion-RecursiveCalibration,
	Fusion12_Fischer,
	Fusion12_Gilitschenski,
	Fusion12_Hanebeck,
	Fusion12_Hekler,
	Fusion12_Noack,
	Fusion12_Packi,
	Fusion12_Reinhardt,
	Fusion12_Reinhardt-FastCI,
	ACC12_Hekler,
	ACC12_Lyons,
	SPIE12_FischerBaum,
	TRB12_Kretz,
	CDC11_Baum,
	CDC11_Hekler,
	ICCV11_Kretz,
	Humanoids11_Rybok,
	SDF11_Huber,
	IROS11_Baum,
	IROS11_Gehrig,
	IROS11_Perez,
	TGF11_Kretz,
	IFAC11_Huber,
	IFAC11_Noack,
	Fusion11_Baum-RHM,
	Fusion11_Baum-USF,
	Fusion11_Baum,
	Fusion11_Bogatyrenko,
	Fusion11_Krauthausen,
	Fusion11_Noack,
	Fusion11_Reinhardt,
	Fusion11_Ruoff,
	% Up tp here 06/24/12
 	ACC11_Bogatyrenko,
 	ACC11_Huber,
 	ACC11_Krauthausen,
 	ACC11_Noack,
	EWSN11_Schmid,
	CDC10_Hekler,
	IROS10_Beutler,
	IROS10_Bogatyrenko,
	IROS10_Packi,
	IROS10_Perez,
	IROS10_Krauthausen,
	SDF10_Baum,
	KI10_Krauthausen,
	IPIN10_Packi,
	VECIMS10_Baum,
	MFI10_Baum,
	MFI10_BaumRect,	
	MFI10_BaumGheta,
	MFI10_Beutler,
	MFI10_Bogatyrenko,
	MFI10_KrauthausenEberhardt,
	MFI10_Krauthausen,
	MFI10_LyonsHekler,
	MSC10_HeklerLyonsNoack,
	Fusion10_BaumKlumpp,
	Fusion10_BaumNoack,
	Fusion10_Beutler,
	Fusion10_DunauPacki,
	Fusion10_EberhardtKlumpp,
	Fusion10_Klumpp-SGMF,
	Fusion10_Klumpp,
	Fusion10_KrauthausenHuber,
	Fusion10_Noack,
	SUTC10_Lyons,
	ACC10_Eberhardt,
	ACC10_Krauthausen,
	PED10_Perez,
%	NIPS09_Deisenroth,
	ISSPIT09_Baum,
%	Deisenroth2009,
	CDC09_HanebeckHuber,
	IROS09_Bogatyrenko,
	SDF09_Baum,
	CCA09_RauhHanebeck,
	Fusion09_Baum,
	Fusion09_Beutler,
	Fusion09_Hoerst-SGMF,
	Fusion09_Huber,
	Fusion09_Klumpp-FusionND,
	Fusion09_Klumpp-Type2Dens,
	Fusion09_Krauthausen,
	Fusion09_Noack,
%	Fusion09_Schieferdecker,
	MSRL09_Deisenroth,
	ICINCO09_Perez,
	ICML09_Deisenroth,
	CASA09_Perez,
	ICASSP09_Beutler,
	MFI09_Beutler,
	ICASSP09_Beutler,
	MFI08_Hanebeck-LCD,
	MFI08_HuberBailey,
	MFI08_Klumpp-DiracFusion,
	MFI08_Klumpp-DiracTrees,
	MFI08_NaWang,
	MFI08_Sawo-Localization,
	MFI08_Schwarzenberg,
	MFI08_WangSzabo,
	MFI08_WeisselSchreiter,
	IFAC08_Huber,
	IFAC08_SawoBeutler,
	IFAC08_WeisselHuber,
	Fusion08_Huber-PGMR,
	Fusion08_Huber_PLSS,
	Fusion08_KlumppSawo-SGMF,
	Fusion08_Noack-ConvexSets,
	Fusion08_SawoKlumpp-SimParam,
	Fusion08_Wang,
	CDC07_HanebeckSchrempf,
	CDC07_WeisselHuber,
	IROS07_SchrempfAlbrecht,
	IROS07_WeisselHuber,
	SDF07_Huber,
	HCI07_Wang,
	Fusion07_HanselmannSchrempf,
	Fusion07_Huber,
	Fusion07_SawoHuber-ParamIdent,
	ACC07_Huber,
	ACC07_Schrempf-DiracMixt,
	ACC07_WeisselHuber,
	ECC07_DeisenrothWeissel,
	CARS07_Roberts,
	ICINCO07_Schrempf-DiracMixt,
	ICINCO07_WeisselHuber,
	IPSN07_Huber,
	ICRA07_Bader,
	WPNC07_Wang,
	CDC06_Brunn-Fourier,
	CDC06_Schrempf-DiracMixt,
	IROS06_Roessler,
	HCRS06_SchrempfSchmid,
	HCRS06_RoesslerSchrempf,
	CCA06_Sawo-ParameterizedGM,
	MFI06_Beutler,
	MFI06_Brunn-Fourier,
	MFI06_DeisenrothWeissel,
	MFI06_Huber,
	MFI06_Schrempf-CramerMises,
	MFI06_Wang,
	ICINCO06_Roessler,
	ICINCO06_SawoRoberts,
	Fusion06_Klumpp,
	Fusion06_Sawo-ParameterizedGMs,
	Fusion06_SchrempfHanselmann,
	ECC05_Rauh,
	INTETAIN05_Groenda,
	CSIT05_Weissel,
	ICINCO05_Hahn,
	ICINCO05_Roessler,
	ICINCO05_Schrempf-Intention,
	ROMAN05_SchrempfSchmid,
%	CCA05_Sawo,
	IROS05_Roessler,
	Fusion05_Roberts,
	Fusion05_Schrempf-MixtureReduction,
	IFAC05_Schrempf-HybridBN,
	ICRA05_Brunn-Splines,
	ICASSP05_Beutler,
	Sensors04_Beutler,
	CSIT04-Schrempf_HybridBN,
	SMC04_Roessler,
	MECHROB04_Roessler,
%	IPDPS04_Weissel,
	CDC03_Feiermann-ProgBayes,
	ECC03_Hanebeck,
	ECC03_Horn-SetBasedLocalization,
	SYSID03_Hanebeck,
	MFI03_Hanebeck,
	SPIE03_HanebeckBriechle-ProgBayes,
	SPIE03_Horn,
	SPIE03_Nitzsche,
	SPIE03_Rauh,
	HCRS02_Nitzsche,
%	DSP02_Roessler,
	ICRA02_Lorch,
	LOCELLUS02_Horn,
	ROMAN01_Nitzsche,
	ECC01_Hanebeck-SetTheoreticEst,
	ECC01_HanebeckHorn-BoundsRandomNoise,
	ECC01_HanebeckHorn-MixStochSetUncertainty,
	MFI01_Hanebeck-PseudoEllipsoids,
	MFI01_HanebeckBriechle-UnknownCorr,
	MFI01_HanebeckHorn-CovarianceBound,
	ACC01_Hanebeck,
	ACC01_HanebeckHorn-MixStochSetUncertainty,
	GMA01_Hanebeck,
	SPIE01_BriechleHanebeck_CrossCorr,
	SPIE01_HanebeckHorn-SLAM,
	RAAD00_Saldic,
	CDC99_Hanebeck,
	SMC99_Briechle-SLAM,
	SMC99_Hanebeck-MixStochSetUncert,
	IROS99_HanebeckSaldic-WheelSystem,
	IROS99_HanebeckHorn-StochSetUncert,
	ECC99_Hanebeck,
	ICRA99_Hanebeck,
	SPIE99_Hanebeck,
	ROMAN98_Schmidt,
	AMC98_Furtwaengler,
	ICRA98_Hanebeck,
	IARP97_Hanebeck,
	IROS97_Hanebeck,
	ISER97_Hanebeck,
	ICASSP96_Hanebeck,
	ICRA96_Hanebeck-SetBayes,
	ICRA96_HanebeckSchmidt-Angle,
	ICASAV95_Hanebeck,
	TELEMAN95_Maier,
	ICONIP94_Hanebeck,
	IROS94_Hanebeck,
	EUFIT94_Hanebeck,
	IFAC93_Hanebeck}
  \putbib
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Artikel in B\"uchern}\English{Papers in Books} (18)}
\begin{bibunit}
  \nocite{
	VMS10_Lyons,
	VMS10_Noack,
	VMS10_Kuwertz,
	LNEE_ICINCO08_WeisselHuber,
	LNEE_ICINCO08_Schrempf,
	LNEE_ICINCO07_SawoRoberts,
	IMS06_BrunnSawo,
	AMS05_RobertsSzabo,
	AMS03_Roessler,
	AMS00_Briechle,
	Robotik00_Hanebeck,
	AMS99_Briechle,
	AMS98_Furtwaengler,
	ISER98_Hanebeck,
	AMS96_Daxwanger,
	AMS95_Hanebeck,
	IRS95_Hanebeck,
	KI95_Hanebeck}
  \putbib
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Artikel in Zeitschriften}\English{Journal Papers} (20)}
\begin{bibunit}
  \nocite{
   	% From here 11/02/12
 	TAES_Baum,
	TRAC12_Deisenroth,  	
  	% SYSAES_Baum,
	% To here 11/02/12
	TM10_Noack,
	TM10_Kuwertz,
	IJCARS10_Bogatyrenko,
	ITM07_BrunnSawo,
	AT05_Briechle,
	AT05_Grigoras,
	PRE04_Nitzsche,
	RA04_Briechle,
	SPL03_Rauh,
	AT03_Nitzsche,
	DGON02_Horn,
	IF00_Hanebeck,
	AT00_Hanebeck,
	AUTO99_Hanebeck,
	JIRSTA98_Ettelt,
	CEP96_Hanebeck,
	FSS96_Hanebeck,
	AT93_Hanebeck}
  \putbib
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Sonderhefte, B\"ucher}\English{Special Issues, Books} (3)}
\begin{bibunit}
  \nocite{
    KI10_Hanebeck,
	AT05_Hanebeck,
	AMS99_Schmidt}
  \putbib
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Habilitation, Dissertation, Diplomarbeit, Studienarbeit}\English{Habilitation, Dissertation, Diploma Thesis, Student Research Project} (4)}
\begin{bibunit}
  \nocite{
	Habil02_Hanebeck,
	Diss97_Hanebeck,
	DA91_Hanebeck,
	SA89_Hanebeck}
  \putbib
\end{bibunit}
ISAS_SADA_Keys.tex
\renewcommand{\refname}{\normalsize \Deutsch{Diplomarbeiten in Karlsruhe}\English{Supervised Diploma Theses in Karlsruhe}}
\begin{bibunit}
  \nocite{
	DA10_Berg,
	DA08_Hirsch,
	DA08_Gehrig,
	DA08_Moellers,
	DA08_Khairallah,
	DA07_Bluemel,
	DA07_Richter,
	DA07_Stiegeler,
	DA07_Schreiter,
	DA07_Krekeler,
	DA07_Meyer,
	DA07_Groene,
	DA07_Seth,
	DA07_Herdtweck,
	DA07_Hanselmann,
	DA07_Albrecht,
	DA07_Krieg,
	DA06_Kern,
	DA06_Beisler,
	DA06_Deisenroth,
	DA06_Klumpp,
	DA06_Huber,
	DA05_Hauth,
	DA05_Prokop,
	DA04_Schrempf,
	DA04_Hahn,
	DA03_Mayer,
	DA03_Misch}
  \putbib%[ISAS_SADA]
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Studienarbeiten in Karlsruhe}\English{Supervised Student Research Projects in Karlsruhe}}
\begin{bibunit}
  \nocite{
	SA10_Reinhardt,
	SA10_Petkov,
	SA08_Kuwertz,
	SA08_Hoerst,
	SA07_Reinbold,
	SA07_Noack,
	SA07_Packi,
	SA07_StefanMeyer,
	SA07_Ota,
	SA07_Zaschke,
	SA07_Stiegeler,
	SA07_Eberhardt,
	SA07_Hauth,
	SA06_Schieferdecker,
	SA06_Baer,
	SA06_Ek,
	SA06_Eberling,
	SA06_Winkler,
	SA06_Stoll,
	SA06_Mende,
	SA06_Hessel,
	SA06_Omid,
	SA06_Armstrong,
	SA06_Hanselmann,
	SA06_Herb,
	SA05_Juetting,
	SA05_Kern,
	SA05_JohannesMeyer,
% 	SA05_GerhardMeyer,	% existiert nicht: ist SA05_JohannesMeyer
	SA05_Krieg,
	SA05_Schamburek,
	SA04_Deisenroth,
	SA04_Lamprecht}
  \putbib%[ISAS_SADA]
\end{bibunit}

\renewcommand{\refname}{\normalsize \Deutsch{Studien--, Semester-- und Diplomarbeiten an der TU M\"unchen}\English{Supervised Student Research Projects and Diploma Theses at TU M\"unchen}}
\begin{bibunit}
  \nocite{
	DA02_Rauh,
	DA01_Riegel,
	DA01_Dabringhaus,
	DA01_Huber,
	SA00_Lohninger,
	SA00_Lohn,
	DA00_Molitor,
	DA99_Russ,
	SA99_Leclerc,
	DA97_Kroessler,
	SA97_Urli,
	DA96_Suess,
	DA95_Maier,
	DA94_Stegmeier,
	DA94_Richter,
	DA93_Berger,
	DA93_Grabl,
	DA93_Boegholm,
	SA93_Roul,
	SA93_Pelat,
	SA92_Faeustle}
  \putbib%[ISAS_SADA]
\end{bibunit}



Deckblaetter.tex
\documentclass[deutsch]{ISASskript_DB}
   
 
   \usepackage[latin1]{inputenc}
   
   
   
   \what{Materialien zur "Ubung}
   \what{Materialien zur Vorlesung}
   
   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   
   \title{Informationsverarbeitung in Sensornetzwerken}
   \author{Prof. Dr.-Ing. Uwe D. Hanebeck,
 %  \\  Dipl.-Ing. Frederik Beutler, 
   \\  Dipl.-Ing. Florian Wei�el}
   \semester{Wintersemester 2007/2008}
  
   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     
%   \title{Stochastische Informationsverarbeitung}
%	 \author{Prof. Dr.-Ing. Uwe D. Hanebeck,\\Dipl.-Inform. Oliver Schrempf,
%	 \\Dipl.-Inform. Peter Krauthausen
%	 }
%	 \semester{Wintersemester 2006/2007}
  
  
%  \title{Lokalisierung mobiler Agenten}
%\author{Prof. Dr.-Ing. Uwe D. Hanebeck,\\ Dipl.-Inform. Marco Huber}

   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%   \title{Unscharfe Mengen}
%   \author{Prof. Dr.-Ing. Uwe D. Hanebeck,\\
%    Dipl.-Inform. Patrick R"o"sler,\\  Dipl.-Inform. Vesa Klumpp}
 %  \semester{Sommersemester 2008}


\begin{document}
\end{document}

DietrichsMatheMacrosDoc.tex
\documentclass[a4paper,11pt,bibtotocnumbered,DIV12]{scrartcl}

\usepackage{ucs} % unicode
\usepackage[utf8x]{inputenc}
\usepackage{ngerman} % deutsche Sprache
\usepackage{bibgerm}
% Schriften
\usepackage[T1]{fontenc} % Output encoding
\usepackage{textcomp}
%\usepackage{pifont} % Symbole(?)
% Times/Helvetica statt Computer modern:
\iffalse % wegen whizzytex
 \usepackage{mathptmx}
 \usepackage[scaled=.90]{helvet}
 \usepackage{courier}
\fi
\renewcommand{\familydefault}{\sfdefault}

\usepackage{xcolor,graphicx}
\usepackage{epstopdf}
\usepackage{tikz}
\usetikzlibrary{arrows}
\usepackage{subfig}

\usepackage[right]{eurosym}
\usepackage{tabularx}

\usepackage{mdwlist}
\usepackage{hyperref}
\usepackage{xspace}
\RequirePackage[squaren,thinspace]{SIunits}

%\usepackage{../ISASPackages/ISASMatheMacros}
\usepackage{../ISASPackages/DietrichsMatheMacros}

\pdfstringdefDisableCommands{\def\\{}} % supresses hyperref warning
\hypersetup{pdftitle={Dokumentation zu der am ISAS verwendeten LaTeX Mathe-Macros},
  pdfsubject={LaTeX Dokumentation},
  pdfauthor={Dietrich Brunn},
  pdfkeywords={LaTeX, Macros, Mathe, Mathematik}, pdfstartview=FitH,
  pdfpagemode={UseOutlines},% None, FullScreen, UseOutlines
  bookmarksnumbered=true, bookmarksopen=false, colorlinks=false
   % a4paper=true,
 }

\title{Dokumentation zu der am ISAS verwendeten \LaTeX Mathe-Macros}
\author{Dietrich Brunn}



\begin{document}
\maketitle
\begin{abstract}
  In diesem Dokument werden die Macros erklärt, welche in den Dateien
  \texttt{ISASMatheMacros.sty} und \texttt{DietrichsMatheMacros.sty} definiert
  sind. Sie enthalten Befehle, welche der am ISAS gebräuchlichen Konventionen für
  Formeln und Formelzeichen entsprechen.
\end{abstract}

\section{Mathematische Macros}
\label{sec:mathematische-macros}


\subsection{Das Aussehen von Variablen betreffend}
\label{sec:das-aussehen-von}

\noindent\begin{tabularx}{1.0\linewidth}{|l|l|X|}
\hline
\LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline
%Variablen
$a$ & \verb+a+ & 
Eine Variable $a$.\\
$\vec{a}$ & \verb+\vec+\verb+{a}+ &
Eine Vektor $\vec{a}$ (Kleinbuchstaben unterstrichen).\\
$\mat{A}$ & \verb+\mat+\verb+{A}+ &
Eine Matrix $\mat{A}$ (fette Großbuchstaben).\\
$\Menge{A}$ & \verb+\Menge+\verb+{A}+ &
Eine Menge $\Menge{A}$ (Großbuchstaben, kalligraphisch).
\\\hline
$\rv{a}$ & \verb+\rv+\verb+{a}+ &
Eine Zufallsvariable $\rv{a}$ (fette Kleinbuchstaben).\\
$\mean{a}$ & \verb+\mean+\verb+{a}+ &
Der Erwartungswert hat einen Hut.\\
$\real{a}$ &\verb+\real+\verb+{a}+ &
Eine Realisierung besitzt einen Stern.%
\\\hline
$\rvec{a}$ & \verb+\rvec+\verb+{a}+ &
Ein Zufallsvektor $\rvec{a}$ (fette Kleinbuchstaben, unterstrichen).\\
$\mean{\vec{a}}$ & \verb+\mean+\verb+{\vec{a}}+ &
Der Erwartungswert hat einen Hut.\\
$\real{\vec{a}}$ &\verb+\real+\verb+{\vec{a}}+ &
Eine Realisierung besitzt einen Stern.\\
\hline
\end{tabularx}

\subsection{Spezielle Variablen und Mengen}
\label{sec:spezielle-mengen}

\noindent\begin{tabularx}{1.0\linewidth}{|l|l|X|}
  \hline \LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline $\vOO, \vI$ &
  \verb+\vOO, \vI+ &
  Der Nullvektor $\vOO= [0\ 0\ \ldots]\T$ und der  Eins-Vektor $\vI= [1\ 1\ \ldots]\T$.\\
  $\mI, \mOO$ & \verb+\mI, \mOO+ &
  Die Einheitsmatrix und die Nullmatrix.
  \\\hline
  $\jj$ & \verb+\jj+ & Die komplexe Einheit ($\jj^2=-1$).  
  \\\hline
  $\IC$ & \verb+\IC+ & Menge der komplexen Zahlen.\\
  $\IN$ & \verb+\IN+ & Menge der natürlichen Zahlen.\\
  $\IR$ & \verb+\IR+ & Menge der reellen Zahlen.\\
  $\IZ$ & \verb+\IZ+ & Menge der ganzen Zahlen.\\
  $\MOO$& \verb+\MOO+ & Die lehre Menge.
  \\\hline
  $\IRp{n}$ & \verb+\IRp+\verb+{n}+ & 
  Ein $n$-dimensionaler Vektorraum.\\
  $\IRp[n]{m}$ & \verb+\IRp+\verb+[n]{m}+ & 
  Ein $n$ mal $m$ dimensionaler Vektorraum.\\
  $\mathds{A}$ & \verb+\mathds+\verb+{A}+ & Ein Buchstabe im
  ,,Blackboard''-Font.  \\\hline
  $\PotMenge{\MA}$& \verb+\PotMenge{\MA}+ & Die Potenzmenge von $\MA$.\\\hline
\end{tabularx}

\subsection{Vordefinierte Variablen}
\label{sec:vord-vari}

\noindent\begin{tabularx}{1.0\linewidth}{|l|l|X|}
\hline
\LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline
% Vektoren:
$\vx, \valpha$ & \verb+\vx, \valpha+& \textbf{Vektoren} sind unterstrichene
Kleinbuchstaben. Sie sind immer stehend, also $\in\IRp[n]{1}$. Es existieren die
Makros $\va$ -- $\vz$ (\verb+\va+ -- \verb+\vz+) und $\valpha$ -- $\vomega$
(\verb+\valpha+ -- \verb+\vomega+).\\
$\vx\T,(\va+\vb)\T$ &
\verb?\vx?\verb?\T,(\va+\vb)\T?&
Transponierte Vektoren sind immer liegend, also $\in\IRp[1]{n}$.\\
$\vOO, \vI$ &
\verb+\vOO, \vI+ &
Der Nullvektor $\vOO= [0\ 0\ \ldots]\T$ und der  Eins-Vektor $\vI= [1\ 1\ \ldots]\T$.
\\\hline
% Matrizen
$\mA,\mLambda$ & \verb+mA, \mLambda+ & \textbf{Matrizen} sind fette
Großbuchstaben. Es existieren die Makros $\mA$ -- $\mZ$ (\verb+\mA+ --
\verb+\MZ+) und $\mGamma$ (\verb+\mGamma+), $\mDelta$ (\verb+\mDelta+), $\mTheta$
(\verb+\mTheta+), $\mLambda$ (\verb+\mLambda+), $\mPsi$ (\verb+\mPsi+), $\mXi$
(\verb+\mXi+), $\mUpsilon$ (\verb+\mUpsilon+), $\mPhi$ (\verb+\mPhi+), $\mOmega$
(\verb+\mOmega+).\\
$\mI, \mOO$ & \verb+\mI, \mOO+ &
Die Einheitsmatrix und die Nullmatrix.\\
$\inv{\mA}, \mA\T$ & \verb+\inv+\verb+{\mA}, \mA\T+ &
Inverse und Transponierte zu  $\mA$.\\
$\adj{A}$ & \verb+\adj+\verb+{A}+& Die Adjungierte zu der Matrix $\mA$.  
\\\hline
%Mengen
$\MA,\ldots,\MZ$ &\verb+\MA,\ldots,\MZ+& \textbf{Mengen} sind kaligraphische
Großbuchstaben. Es exisiert $\MA$ bis $\MZ$.  \\\hline
% Zufallsvariablen:
$\rra, \rrxi$ &
\verb+\rra, \rrxi+ &
\textbf{Zufallsvariablen} sind fette Kleinbuchstaben.\\
$\eeb, \eegamma$  & 
\verb+\eeb, \eegamma+&
Erwartungswerte besitzen einen Hut.\\
$\rlb, \rlgamma$  & 
\verb+\rlb, \rlgamma+&
Realisierungen besitzen einen Stern.\\
\hline
% Zufallsvektoren:
$\rva, \rvxi$ &
\verb+\rva, \rvxi+ &
\textbf{Zufallsvektoren} sind unterstrichene fette Kleinbuchstaben.\\
$\evb, \evgamma$  & 
\verb+\evb, \evgamma+&
Erwartungswerte besitzen einen Hut.\\
$\rlvb, \rlvgamma$  & 
\verb+\rlvb, \rlvgamma+&
Realisierungen besitzen einen Stern.
\\\hline
\end{tabularx}

\subsection{Spezielle Funktionen}
\label{sec:spezielle-funktionen}

\noindent\begin{tabularx}{1.0\linewidth}{|l|l|X|}
\hline
\LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline
$\Pr(\rrx=x)$ & \verb+\Pr+\verb+(\rrx+\verb+=x)+ & Die Wahrscheinlichkeit, dass $\rrx=x$ ist.\\
$\fd{\rrx}(x)$ & \verb+\fd+\verb+{\rrx}(x)+ &
Die Wahrscheinlichkeitsdichte der Zufallsvariable $\rrx$.\\
$\Fd{\rrx}(x)$ & \verb+\Fd+\verb+{\rrx}(x)+ &
Die kumulative Wahrscheinlichkeitsdichte $\Pr(\rrx\leq x)$ der Zufallsvariable $\rrx$.
\\
$\Pd{\rrx}(x)$ & \verb+\Pd+\verb+{\rrx}(x)+ &
Die $\Psi$-Dichte der Zufallsvariable $\rrx$.\\\hline
$\fd[h]{\rvx}(\vx)$ & \verb+\fd+\verb+[h]{\rvx}(\vx)+ &
Die Wahrscheinlichkeitsdichte der Art h des Zufallsvektors $\rrx$.\\
$\Fd[h]{\rvx}(\vx)$ & \verb+\Fd+\verb+[h]+\verb+{\rvx}(\vx)+ &
Die kumulative Wahrscheinlichkeitsdichte $\Pr(\rvx\leq \vx)$ der Art h des Zufallsvektors $\rrx$.
\\
$\Pd[h]{\rvx}(\vx)$ & \verb+\Pd+\verb+[h]+\verb+{\rvx}(\vx)+ &
Die $\Psi$-Dichte der Art h des Zufallsvektors $\rrx$.
\\\hline
$\fde{\rrx}(x)$ & \verb+\fde+\verb+{\rrx}(x)+ &
Die geschätzte (estimated) Wahrscheinlichkeitsdichte der Zufallsvariable $\rrx$.\\
$\Fde{\rrx}(x)$ & \verb+\Fde+\verb+{\rrx}(x)+ &
Die geschätzte (estimated) kumulative Wahrscheinlichkeitsdichte $\Pr(\rvx\leq
\vx)$ der Zufallsvariable $\rrx$.
\\
$\Pde{\rrx}(x)$ & \verb+\Pde+\verb+{\rrx}(x)+ &
Die geschätzte (estimated) $\Psi$-Dichte der Zufallsvariable $\rrx$.
\\\hline
$\fdL{\rrx}(x)$ & \verb+\fde+\verb+{\rrx}(x)+ &
Die Likelihood der Zufallsvariable $\rrx$.\\
$\FdL{\rrx}(x)$ & \verb+\Fde+\verb+{\rrx}(x)+ &
Die geschätzte kumulative Likelihood $\Pr(\rvx\leq \vx)$ der Zufallsvariable $\rrx$.
\\
$\PdL{\rrx}(x)$ & \verb+\Pde+\verb+{\rrx}(x)+ &
Die geschätzte $\Psi$-Dichte der Zufallsvariable $\rrx$.
\\\hline
$\fdp{\rrx}(x)$ & \verb+\fdp+\verb+{\rrx}(x)+ &
Die posteriore Wahrscheinlichkeitsdichte der Zufallsvariable $\rrx$.\\
$\Fdp{\rrx}(x)$ & \verb+\Fdp+\verb+{\rrx}(x)+ &
Die posteriore kumulative Wahrscheinlichkeitsdichte $\Pr(\rvx\leq \vx)$ der Zufallsvariable $\rrx$.
\\
$\Pdp{\rrx}(x)$ & \verb+\Pdp+\verb+{\rrx}(x)+ &
Die posteriore $\Psi$-Dichte der Zufallsvariable $\rrx$.
\\\hline
$\fdT{\rrx}(x)$ & \verb+\fdT+\verb+{\rrx}(x)+ &
Die Transitionsdichte der Zufallsvariable $\rrx$.\\
$\FdT{\rrx}(x)$ & \verb+\FdT+\verb+{\rrx}(x)+ &
Die kumulative Transitionsdichte $\Pr(\rvx\leq \vx)$ der Zufallsvariable $\rrx$.
\\
$\PdT{\rrx}(x)$ & \verb+\PdT+\verb+{\rrx}(x)+ &
Die  Transitions-$\Psi$-Dichte der Zufallsvariable $\rrx$.
\\\hline
\end{tabularx}


\subsection{Operatoren (ohne Argumente)}
\label{sec:oper-ohne-argum}

\noindent\begin{tabularx}{1.0\linewidth}{|l|l|X|}
\hline
\LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline
$\dd$ & \verb+\dd+ & 
Der Differential-Operator, wie in $\dd/\dd x$ oder in $\int\dd \xi$.\\
$\partial$ &\verb+\partial+ &
Der partielle Differential-Operator.\\
$(.)\T$ & \verb+(.)\T+ &
Der Transponiert-Operator.
\\\hline
$\argmin$ & \verb+\argmin+ &
Argument des Minimums.\\
$\argmax$ & \verb+\argmax+ &
Argument des Maximums.\\
$\atan2$ & \verb+\atan+\verb+2+ & Arcus-Tangens mit zwei Argumenten.\\
$\diag$ & \verb+\diag+ &
 Diagonale einer Matrix\\
$\erf$ & \verb+\erf+ & Fehler-Funktion.\\
$\Gauss$ &\verb+\Gauss+ & Gaußdichte.\\
$\rang$ & \verb+\rang+ & Rang einer Matrix (deutsch).\\
$\rank$ & \verb+\rank+ & Rang einer Matrix (englisch).\\
\hline
$\rect_b(x)$ & \verb+\rect+\verb+_b(x)+ & Rechteckfunktion:
 \begin{equation*}
   \rect_b(x) :=
  \begin{cases}
    1   & -\tfrac{b}{2} < x < \tfrac{b}{2} \\
    1/2 & |x| = \tfrac{b}{2}\\
    0   & \text{sonst}
  \end{cases}
 \end{equation*}\\
 $\step(x)$ & \verb+\step(x)+ & Einheitssprung:
\begin{equation}
 \step(x) :=
 \begin{cases}
   0   & x < 0\\
   1/2 & x = 0\\
   1   & x > 0
 \end{cases}
\end{equation}\\
$\sign(x)$ & \verb+\sign(x)+& Signumfunktion.\\\hline
 $\Ind{\MA}(x)$ &\verb+\Ind+\verb+{\MA}(x)+& Die Indikatorfunktion
\begin{equation}
 \Ind{\MA}(x) :=  
 \begin{cases}
    1 & x\in\MA\ ,\\
    0 & x\notin\MA\ .
   \end{cases}
\end{equation}
\\\hline
\end{tabularx}

\subsection{Operatoren (mit Argumente)}
\label{sec:oper-ohne-argum}

\noindent\begin{tabularx}{1.0\linewidth}{|l|l|X|}
\hline
\LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline
%
$\inv{.}$ & \verb+\inv+\verb+{.}+ &
Der Inversions-Operator.\\
$\adj{(.)}$ & \verb+\adj+\verb+{(.)}+ &
Der Adjungiert-Operator
\\\hline
$\Re[z]$ & \verb+\Re+\verb+[z]+ & Realteil einer Zahl.\\
$\Im[z]$ & \verb+\Im+\verb+[z]+ & Imaginärteil einer Zahl.\\
$\E{\rrx}$ & \verb+\E+\verb+{\rrx}+ & 
Erwartungswert einer Zufallsvariablen.
\\\hline%
$\E[\fd{x}]{\rrx}$ & \verb+\E+\verb+[\fd+\verb+{x}]{\rrx}+ & 
Erwartungswert einer Zufallsvariablen mit der Dichte $\fd{x}(x)$.\\
$\Cov[\fd{\rvx}]{\rvx}$ & \verb+\Cov+\verb+[\fd+\verb+{\rvx}]{\rvx}+ &
Die Kovarianz des Zufallsvektors $\rvx\sim\fd{\rvx}(\rvx)$.\\
$\Var[\fd{\rvx}]{\rvx}$ & \verb+\Var+\verb+[\fd{\rvx}]{\rvx}+ &
Die Varianz des Zufallsvektors $\rvx\sim\fd{\rvx}(\rvx)$.
\\\hline
$\diff{x} h(x)$ & \verb+\diff+\verb+{x} h(x)+ &
Die Ableitung von $h(x)$.\\
$\diff[h]{x}$ & \verb+\diff+\verb+[h]{x}+ &
Die Ableitung von $h$ nach $x$.\\
$\pardiff{x} h(x)$ & \verb+\pardiff+\verb+{x} h(x)+ &
Die partielle Ableitung von $h(x)$.\\
$\pardiff[h]{x}$ & \verb+\diff+\verb+[h]{x}+ &
Die partielle Ableitung von $h$ nach $x$.
\\\hline
$\sprod{f}{g}$& \verb+\sprod+\verb+{f}{g}+ & Das Skalarprodukt der Funktionen $f$
und $g$
\\\hline
$\floor{x}$& \verb+\floor+\verb+{x}+ & Abrunden von $x$ auf die nächste ganze
Zahl.\\
$\ceil{x}$& \verb+\ceil+\verb+{x}+ & Aufrunden von $x$ auf die nächste ganze Zahl.
\\\hline
\end{tabularx}


\subsection{Setzen von Einheiten (\texttt{SIunits})}
\label{sec:setzen-von-einheiten}

Mit dem Paket \texttt{SIunits}  können SI-Einheiten
normgerecht gesetzt werden. Empfohlen wird
%\begin{verbatim}\usepackage[squaren,thinspace]{SIunits}\end{verbatim}
zu verwenden. Ferner wird sicher gestellt, dass der Abstand zwischen Zahl und
Einheit korrekt ist. Weitere Infos befinden sich in der zugehörigen Doku vom
Paket \texttt{SIunits}.\par

\noindent\begin{tabularx}{1.0\linewidth}{|l|X|X|}
\hline
\LaTeX & \LaTeX-Code & Beschreibung\\\hline\hline
%
$\unit{1}\centi\meter$ & 
\verb+\unit+\verb+{1}\centi+\verb+\meter+
& Einheiten werden mit
\verb|\unit| abgesetzt.\\\hline
\unit{1}\kilogram, \unit{1}\ampere, \unit{1}\watt\ldots &
\verb+\unit+\verb+{1}\kilogram+, \verb+\unit+\verb+{1}\ampere+, \verb+\unit+\verb+{1}\watt+ \ldots & 
Fast alle Einheiten sind im  Paket \texttt{SIunits} vordefiniert (siehe
zugehörige Doku)\\\hline
\unit{1}\nano\meter, \unit{1}\micro\meter, \ldots &
\verb+\unit+\verb+{1}\nano+\verb+\meter+, \verb+\unit{1}+\verb+\micro+\verb+\meter+, \ldots &
die gängingen SI-Präfixe ebenfalls.\\\hline
{\volt\ampere}, {\meter\per\second}&
\verb+{\volt\ampere}+, \verb+{\meter\per\second}+&
Es können  auch Einheiten kombiniert werden\\\hline
{\centi\cubic\meter}, {\rpsquare\second} & 
\verb+{\cubic\centi\meter}+, \verb+{\rpsquare\second}+
& Es gibt auch quadrierte und kubische Einheiten\\\hline
$\color{red}x\ [\ohm]$ && \textbf{\color{red}Falsch!!} [.] bedeutet \emph{Einheit von}\\\hline
$[x] = \ohm$&\verb+[x] = \ohm+ & Richtig\\\hline
\end{tabularx}
\end{document}

%%% Local Variables: 
%%% mode: latex
%%% TeX-master: t
%%% End: 

% LocalWords:  kalligraphisch Blackboard Inversions

appendix.tex
\section{Proof of Lemma \ref{lemma:partialderivatives}.}
\label{sec:proofpartialderivatives}
\noindent
%\begin{proof}
We use the identities
\[ {}_1F_1 \left(\frac{1}{2}, 1, \begin{pmatrix} z_1 & 0 \\ 0 & z_2 \end{pmatrix} \right) =  \exp(z_2) \cdot {}_1F_1 \left( \frac{1}{2},1,z_1-z_2 \right) \]
and
\[ \frac{\partial}{\partial z} {}_1F_1 (a,b,z) = \frac{a}{b} {}_1F_1 (a+1, b+1, z) \ . \]
This yields
\begin{align*}
\frac{\partial}{\partial z_1} \frac{F}{|S_{d-1}|}
& = \frac{\partial}{\partial z_1} {}_1F_1 \left( \frac{1}{2}, 1, \begin{pmatrix} z_1 & 0 \\ 0 & z_2 \end{pmatrix}\right)\\ 
& = \exp(z_2) \frac{\partial}{\partial z_1} {}_1F_1 \left( \frac{1}{2}, 1, z1-z2 \right) \\ 
& = \exp(z_2) \frac{1}{2} {}_1F_1 \left( \frac{3}{2}, 2, z1-z2 \right)
\end{align*}
and
\begin{align*}
\frac{\partial}{\partial z_2} \frac{F}{|S_{d-1}|}
= & \frac{\partial}{\partial z_2} {}_1F_1 \left( \frac{1}{2}, 1, \begin{pmatrix} z_1 & 0 \\ 0 & z_2 \end{pmatrix}\right) \\
= &\frac{\partial}{\partial z_2} \left( \exp(z_2) {}_1F_1 \left( \frac{1}{2}, 1, z_1 - z_2 \right) \right) \\
= &\exp(z_2) {}_1F_1 \left( \frac{1}{2}, 1, z_1 - z_2 \right) + \exp(z_2) \frac{\partial}{\partial z_2} {}_1F_1 \left( \frac{1}{2}, 1, z_1 - z_2 \right) \\
= &\exp(z_2) {}_1F_1 \left( \frac{1}{2}, 1, z_1 - z_2 \right) - \exp(z_2) \frac{1}{2} {}_1F_1 \left( \frac{3}{2}, 2, z_1 - z_2 \right) \\
= &\exp(z_2) \Bigg( {}_1F_1 \left( \frac{1}{2}, 1, z_1 - z_2 \right) - \frac{1}{2} {}_1F_1 \left( \frac{3}{2}, 2, z_1 - z_2 \right) \Bigg)
\end{align*}
\hfill\qedsymbol
%\end{proof}

\section{Proof of Lemma \ref{lemma:composition}.}
\label{sec:proofcomposition}
\noindent

%\begin{proof}
The covariance of the composition
\begin{align*}
\mat{C} =& \Cov (\vec{x} \oplus \vec{y}) \\
=& \Cov \left( \begin{pmatrix} x_1y_1 - x_2 y_2\\ x_1 y_2 + x_2 y_1 \end{pmatrix} \right) \\
=& \begin{pmatrix}
\Var (x_1y_1 - x_2 y_2) & \Cov(x_1y_1 - x_2 y_2,  x_1 y_2 + x_2 y_1)\\
* & \Var( x_1 y_2 + x_2 y_1)
\end{pmatrix}
\end{align*}
can be obtained by calculating the matrix entries individually. For the first entry we get
\begin{align}
c_{11} =& \Var (x_1y_1 - x_2 y_2) \notag \\
=& \E ((x_1y_1 - x_2 y_2)^2) - (\E(x_1y_1 - x_2 y_2))^2 \notag \\
=& E(x_1^2y_1^2 - 2 x_1y_1 x_2y_2 + x_2^2y_2^2) %\notag \\
- (\E(x_1y_1) - \E(x_2y_2))^2 \label{eq:independence} \\
=& \E(x_1^2)\E(y_1^2) - 2 \E(x_1x_2)\E(y_1y_2) + \E(x_2^2)\E(y_2^2) \label{eq:linearity} \\
&- (\underbrace{\E(x_1)}_{0} \underbrace{\E(y_1)}_{0} - \underbrace{\E(x_2)}_{0} \underbrace{\E(y_2)}_{0})^2 \label{eq:zero} \\
=& a_{11} b_{11} - 2 a_{12} b_{12} + a_{22} b_{22} \notag .
\end{align}
We use independence of $\vec{x}$ and $\vec{y}$ in (\ref{eq:independence}), linearity of the expectation value in (\ref{eq:linearity}) and symmetry of the Bingham in (\ref{eq:zero}).
Analogously we calculate
\begin{align*}
c_{22} =& a_{11} b_{22} + 2 a_{12} b_{12} + a_{22} b_{11} \ .
\end{align*}
The off-diagonal entry can be calculated similarly
\begin{align*}
c_{12} =& \Cov(x_1y_1 - x_2 y_2,  x_1 y_2 + x_2 y_1) \\
=& \E((x_1y_1 - x_2 y_2) \cdot (x_1 y_2 + x_2 y_1)) 
- \E(x_1y_1 - x_2 y_2)\cdot \E(x_1 y_2 + x_2 y_1)\\
=& \E (x_1^2y_1y_2 - x_1x_2y_2^2 + x_1x_2y_1^2 - x_2^2y_1y_2 )
- (\E(x_1)\E(y_1) - \E(x_2) \E(y_2)) \\ 
&\cdot (\E(x_1) \E( y_2) + \E(x_2) \E(y_1)) \\
=& a_{11} b_{12} - a_{12} b_{22} + a_{12} b_{11} - a_{22} b_{12} \ .
\end{align*}
\hfill\qedsymbol
%\end{proof}


bingham_distribution.tex
The Bingham distribution appears naturally when a $d$-dimensional normal random vector ${\bf \vec{x}}$ with $\E({\bf \vec{x}})=\vec{0}$ is conditioned on $||{\bf \vec{x}}||=1$ \cite{kume2005}. In the following, we will introduce the Bingham distribution and derive the formulas for multiplication of two Bingham probability density functions. Furthermore, we will present a method for computing the composition of two Bingham-distributed random variables, which is analogous to the addition of real random variables.

% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
\subsection{Probability Density Function}
% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
\begin{definition}
Let $S_{d-1} = \{ \vec{x} \in \mathbb{R}^d : ||\vec{x}||=1 \} \subset \mathbb{R}^d$ be the unit hypersphere in $\mathbb{R}^d$. The probability density function (pdf) 
\[ f: S_{d-1} \to \mathbb{R} \]
of a Bingham  distribution \cite{bingham1974} is given by
\[ f(\vec{x}) = \frac{1}{F} \cdot \exp ( \vec{x}^T \mat{M}\, \mat{Z}\, \mat{M}^T \vec{x} )\ ,  \]
where $\mat{M} \in \mathbb{R}^{d \times d}$ is an orthogonal matrix ($\mat{M}\,\mat{M}^T = \mat{M}^T\,\mat{M} = \mat{I}_{d \times d}$) describing the orientation, $\mat{Z} = \text{diag} (z_1, \dots z_{d-1}, 0) \in \mathbb{R}^{d \times d}$ with $z_1 \leq \dots \leq z_{d-1} \leq 0$ is the concentration matrix, and $F$ is a normalization constant. 
\end{definition}
As Bingham showed, adding a multiple of the identity matrix $ \mat{I}_{d \times d}$ to $\mat{Z}$ does not change the distribution. Thus, we conveniently force the last entry of \mat{Z} to be zero. Because it is possible to swap columns of $\mat{M}$ and the according diagonal entries in $\mat{Z}$ without changing the distribution, we can enforce $z_1 \leq \dots \leq z_{d-1}$. This representation allows us to obtain the mode of the distribution very easily by taking the last column of $\mat{M}$.

The pdf is antipodally symmetric, i.\,e., $f(\vec{x}) = f(- \vec{x})$ holds for all $x \in S_{d-1}$. Consequently, the Bingham distribution is invariant to rotations by $180^\circ$. Examples of the pdf for two dimensions ($d=2$) are shown in Fig. \ref{fig:pdf3d} and Fig. \ref{fig:pdf2d}. The Bingham distribution is very similar to a Gaussian if and only if the uncertainty is small. This can be seen in Fig. \ref{fig:kld}, which shows the Kullback-Leibler divergence between a Bingham pdf and a corresponding Gaussian pdf.  

\begin{figure}[t]
\centering
\includegraphics[width=0.6\textwidth]{figures/binghamPdf3d.eps}
\caption{Bingham pdf with  $\mat{M}=\mat{I}_{2 \times 2}$ and $\mat{Z} = \diag(-8, 0)$ as a 3D plot. This corresponds to a standard deviation of $16^\circ$.}
\label{fig:pdf3d}
\end{figure}

\begin{figure}[t]
\centering
\includegraphics[width=0.6\textwidth]{figures/binghamPdf2d.eps}
\caption{Bingham pdf with $\mat{M}=\mat{I}_{2 \times 2}$ for different values of $\mat{Z} = \diag (z_1, 0)$ and $x=(\cos(\theta), \sin(\theta))^T$. These values for $z_1$ correspond to standard deviations of approximately $36^\circ$, $16^\circ$, and $6^\circ$ respectively.}
\label{fig:pdf2d}
\end{figure}

\begin{figure}[t]
\centering
\includegraphics[width=0.6\textwidth]{figures/kld.eps}
\caption{Kullback-Leibler divergence on the interval $[0, \pi]$ between a Bingham pdf with $\mat{M}=\mat{I}_{2 \times 2}$, $\mat{Z} = \diag (z_1, 0)$ and a Gaussian pdf with equal mode and standard deviation. For small uncertainties ($z_1 < - 15$, which corresponds to a standard deviation of about $11^\circ$), the Gaussian and Bingham distributions are almost indistinguishable. However, for large uncertainties, the Gaussian approximation becomes quite poor.}
\label{fig:kld}
\end{figure}

%todo calculation of bounds?

% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
\subsection{Normalization Constant}
% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
The normalization constant can be calculated with the help of the hypergeometric function of a matrix argument \cite{herz1955,koev2006,muirhead1982}. It is given by
\[ F := |S_{d-1}| \cdot {}_1F_1 \left( \frac{1}{2}, \frac{d}{2}, \mat{Z} \right)\ , \]
where $|S_{d-1}|$ is the surface area of the $d$-sphere and ${}_1F_1(\cdot,\cdot,\cdot)$ is the hypergeometric function of matrix argument. In the two-dimensional case ($d=2$), this reduces to
\[ F = 2 \pi \cdot {}_1F_1 \left( \frac{1}{2}, 1, \begin{pmatrix} z_1 & 0 \\ 0 & 0 \end{pmatrix} \right) = 2 \pi \cdot {}_1F_1 \left(\frac{1}{2}, 1, z_1 \right)\ , \]
so it is sufficient to compute the hypergeometric function of a scalar argument, which is described in \cite{abramowitz1964}.
%todo numerical algorithms for hypergeometric function

% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
\subsection{Multiplication}
For two given Bingham densities, we want to obtain their product. This product is used for Bayesian inference involving Bingham distributions. The result presented below yields a convenient way to calculate the product of Bingham distributions.
\begin{lemma}
Bingham distributions are closed under multiplication with renormalization.
\end{lemma}
\begin{proof}
Consider two Bingham distributions
\[ f_1(\vec{x}) = F_1 \cdot \exp ( \vec{x}^T \mat{M}_1\, \mat{Z}_1\, \mat{M}_1^T \vec{x} )  \]
and 
\[ f_2(\vec{x}) = F_2 \cdot \exp ( \vec{x}^T \mat{M}_2\, \mat{Z}_2\, \mat{M}_2^T \vec{x} ) \ . \] 
Then 
\begin{align*} 
f_1(\vec{x}) \cdot f_2(\vec{x}) 
&= F_1 F_2 \cdot \exp (\vec{x}^T ( \underbrace{\mat{M}_1 \mat{Z}_1 \mat{M}_1^T + \mat{M}_2 \mat{Z}_2 \mat{M}_2^T)}_{ =: \mat{C}} \vec{x}) \\
&\propto F \cdot \exp ( \vec{x}^T \mat{M}\, \mat{Z}\, \mat{M}^T \vec{x} )
\end{align*}
with $F$ as the new normalization constant after renormalization, $\mat{M}$ are the unit eigenvectors of $\mat{C}$, $\mat{D}$ has the eigenvalues of $\mat{C}$ on the diagonal (sorted in ascending order) and $\mat{Z} = \mat{D} - \mat{D}_{dd}\mat{I}_{d \times d}$ where $ \mat{D}_{dd}$ refers to the bottom right entry of $\mat{D}$, i.\,e., the largest eigenvalue.
\end{proof}

% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
\subsection{Estimation of Parameters}
% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
Estimating parameters for the Bingham distribution is not only motivated by the need to estimate noise parameters from samples. It also plays a crucial role in the prediction process when computing the composition of two Bingham random vectors. This procedure is based on matching covariance matrix. Be aware that although the Bingham distribution is only defined on $S_{d-1}$, we can still compute its covariance in $\mathbb{R}^d$. Thus, we will present both the computation of the covariance matrix of a Bingham distributed random vector and the computation of parameters for a Bingham distribution with a given covariance (which could originate from an arbitrary distribution on the hypersphere). 


%\begin{definition}
%Given a matrix of $n$ samples $\mat{X} \in \mathbb{R}^{d \times n}$ with one sample per column, the scatter matrix is given by $\mat{S} = \frac{1}{n} \mat{X} \mat{X}^T \in \mathbb{R}^{d \times d}$. 
%\end{definition}

The maximum likelihood estimate for the parameters $(\mat{M},\mat{Z})$ of a Bingham distribution can be obtained as described in \cite{bingham1974}. $\mat{M}$ can be obtained as the matrix of eigenvectors of the covariance $\mat{S}$ with eigenvalues $\omega_1 \leq \omega_2$. In other words, $\mat{M}$ can be found as the eigendecomposition of $\mat{S} = \mat{M} \cdot \diag (\omega_1, \omega_2) \cdot \mat{M}^T$. To calculate $\mat{Z}$, the equations
\[ \frac{\frac{\partial}{\partial z_i} {}_1F_1 \left( \frac{1}{2}, 1, \begin{pmatrix} z_1 & 0 \\ 0 & z_2 \end{pmatrix} \right)}{{}_1F_1 (\frac{1}{2}, 1, z_1 )}  = \omega_i, \quad i=1,2 \]
have to be solved under the constraint $z_2=0$, which is justified by the argumentation above and used to simplify the computation.  This operation is performed numerically.

Conversely, for a given a Bingham distribution $(\mat{M},\mat{Z})$, the covariance matrix can be calculated according to 
\begin{align*} 
\mat{S} &=  \mat{M} \cdot \diag (\omega_1, \omega_2) \cdot \mat{M}^T \\
&=  \mat{M} \cdot \diag \left( \frac{1}{F} \frac{\partial F}{\partial z_1}, \frac{1}{F} \frac{\partial F}{\partial z_2} \right) \cdot \mat{M}^T  \\
&=\sum_{i=1}^2 \frac{1}{F} \frac{\partial F}{\partial z_i} \mat{M}(:,i) \mat{M}(:,i)^T \ ,
\end{align*}
where $\mat{M}(:,i)$ refers to the $i$-th column of $\mat{M}$ \cite{Glover13}. Thus, any Bingham distribution is uniquely defined by its covariance matrix and vice versa. The following Lemma simplifies the computation of partial derivatives of a confluent hypergeometric function of a $2\times 2$ matrix argument, which is used in computation of the covariance matrix as derived above.

\begin{lemma}
\label{lemma:partialderivatives}
For $d=2$, the partial derivatives
\[  \frac{\partial}{\partial z_i} {}_1F_1 \left( \frac{1}{2}, 1, \begin{pmatrix} z_1 & 0 \\ 0 & z_2 \end{pmatrix}\right) \ , \quad i=1,2 \]
can be reduced to hypergeometric functions of scalar argument.
\end{lemma}
\begin{proof}
See \ref{sec:proofpartialderivatives}.
\end{proof}


% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
\subsection{Composition}
% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
Now, we want to derive the composition of Bingham distributed random variables, which is the directional analogue to adding random variables. This operation can, for example, be used to disturb an uncertain Bingham-distributed system state with Bingham-distributed noise, similar to using a convolution to disturb a probability distribution on $\mathbb{R}$ with additive noise. First, we define a composition of individual points on the hypersphere $S_{d-1}$, which we then use to derive the composition of Bingham distributions.


%\begin{figure}
%\setlength{\tabcolsep}{3pt}
%\begin{tabular}{lccccccc}
%Bingham pdf & $f_\mat{A}$ & $\oplus$ & $f_\mat{B}$ & $=$ & $f_\mat{A} \oplus f_\mat{B}$  & $\approx$ & $f_\mat{C}$ \\
% & $\Downarrow$ & & $\Downarrow$ & & & & $\Uparrow$ \\
%scatter matrix \quad & $\mat{A}$ & $\oplus$ & $\mat{B}$ & $=$ & & & $\mat{C}$ \\
% & $\updownarrow$ & & $\updownarrow$ & & & & $\updownarrow$ \\
%samples & $\mat{X}$ & $\oplus$ & $\mat{Y}$ & $=$ & &  & $\mat{Z}$
%\end{tabular}
%\centering
%\caption{Approximation of composition operation based on scatter matrices.}
%\label{fig:compositionapproximation}
%\end{figure}

The composition of two Bingham distributions depends on the interpretation of the unit vectors, for example as complex numbers or quaternions. We assume that a composition function
\[ \oplus: S_{d-1} \times S_{d-1} \to S_{d-1} \]
is given. The function $\oplus$ has to be compatible with $180^\circ$ degree symmetry, i.e., 
\[ \pm(x \oplus y) = \pm( (-x) \oplus y ) = \pm( x \oplus (-y) ) = \pm( (-x) \oplus (-y) ) \]
for all $x,y \in S_{d-1}$. Furthermore, we require the quotient $(S_{d-1} / \{ \pm 1 \}, \oplus)$ to have an algebraic group structure. This guarantees associativity, the existence of an identity element, and the existence of inverse elements.
%todo group structure, equivalence classes

In the complex case, we interpret $S_1 \subset \mathbb{R}^2$ as unit vectors in $\mathbb{C}$, where the first dimension is the real part and the second dimension the imaginary part. In this interpretation, the Bingham distributions can be understood as a distribution on a subset of the complex plane, namely the unit circle.
\begin{definition}
The composition function $\oplus$ is defined to be complex multiplication, i.e.,
\[ \begin{pmatrix} x_1 \\ x_2\end{pmatrix} \oplus \begin{pmatrix} y_ 1\\ y_2\end{pmatrix} = \begin{pmatrix} x_1y_1 - x_2 y_2\\ x_1 y_2 + x_2 y_1 \end{pmatrix} \]
analogous to
\[ (x_1 + i x_2) \cdot (y_1 + i y_2) = (x_1 y_1 - x_2 y_2) + i (x_1 y_2 + x_2 y_1) \ . \]
\end{definition}
Since we only consider unit vectors, the composition $\oplus$ is equivalent to adding the angles of both complex numbers when they are represented in polar form. The identity element is $\pm 1$ and the inverse element for $(x_1,x_2)$ is the complex conjugate $\pm(x_1,-x_2)$.


% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
Unfortunately, the Bingham distribution is not closed under this kind of composition. That is, the resulting random vector is not Bingham distributed. Thus, we propose a technique to approximate a Bingham distribution to the composed random vector. The composition of two Bingham distributions $f_\mat{A}$ and $f_\mat{B}$ is calculated by considering the composition of their covariance matrices $\mat{A}, \mat{B}$ and estimating the parameters of $f_\mat{C}$ based on the resulting covariance matrix. Composition of covariance matrices can be derived from the composition of random vectors.

\begin{lemma}
Let $f_\mat{A}$ and $f_\mat{B}$ be Bingham distributions with covariance matrices \label{lemma:composition} 
\[\mat{A} = \begin{pmatrix}a_{11} & a_{12} \\ * & a_{22}\end{pmatrix} \text{ and } \mat{B}=\begin{pmatrix}b_{11} & b_{12} \\ * & b_{22}\end{pmatrix}\ , \]
respectively. Let $\vec{x},\vec{y} \in S_1 \subset \mathbb{R}^2$ be independent random vectors distributed according to $f_\mat{A}$ and $f_\mat{B}$. Then the covariance
 \[ \mat{C} =  \begin{pmatrix}c_{11} & c_{12} \\ * & c_{22}\end{pmatrix} :=  \Cov (\vec{x} \oplus \vec{y}) \]
of the composition is given by
\begin{align*}
c_{11} =& a_{11} b_{11} - 2 a_{12} b_{12} + a_{22} b_{22}\ , \\
c_{12} =& a_{11} b_{12} - a_{12} b_{22} + a_{12} b_{11} - a_{22} b_{12}\ , \\
c_{22} =& a_{11} b_{22} + 2 a_{12} b_{12} + a_{22} b_{11}\ .\\
\end{align*}
\end{lemma}
\begin{proof}
See \ref{sec:proofcomposition}.
\end{proof}

Based on $\mat{C}$, maximum likelihood estimation is used to obtain the parameters $\mat{M}$ and $\mat{Z}$ of the uniquely defined Bingham distribution with covariance $\mat{C}$ as described above. This computation can be done in an efficient way, because the solution of the equation involving the hypergeometric function is the only part which is not given in closed form. This does not present a limitation to the proposed algorithm, because there are many efficient ways for the computation of the confluent hypergeometric function of a scalar argument \cite{luke1977, muller2001}.


% % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % %
%Unfortunately, the Bingham distribution is not closed under this kind of composition. That is, the resulting random vector is not Bingham-distributed. Thus we propose a technique to approximate a Bingham distribution to the composed random vector. The basic idea of this method is depicted in Fig. \ref{fig:compositionapproximation}. The composition of two Bingham distributions $f_\mat{A}$ and $f_\mat{B}$ is calculated by considering the composition of their scatter matrices $\mat{A}, \mat{B}$ and estimating the parameters of $f_\mat{C}$ based on the resulting scatter matrix $\mat{C} = \mat{A} \oplus \mat{B}$. Composition of scatter matrices can be derived from the composition of samples. These samples are only needed as part of the derivation and they do not need to be calculated explicitly.
%
%In order to derive the composition of two Bingham distributions we consider the composition of samples corresponding to those distributions. Let $\mat{A}$ and $\mat{B}$ be the scatter matrices of two Bingham distribution. Let us further assume, that $\mat{A}$ and $\mat{B}$ are generated by $n$ samples $\mat{X}$ and $m$ samples $\mat{Y}$, i.\,e. 
%\[ \mat{A} = \frac{1}{n} \mat{X} \mat{X}^T, \mat{B} = \frac{1}{m} \mat{Y} \mat{Y}^T \ , \] even though we do not know or even need these samples explicitly. We want to calculate the scatter matrix 
%\[ \mat{C} = \frac{1}{n\cdot m} \mat{Z} \mat{Z}^T \]
%of the distribution which is estimated from samples $\{\vec{x} \oplus \vec{y} : \vec{x} \in \mat{X}, \vec{y} \in \mat{Y} \}$ which are stored in the matrix $\mat{Z}$. For $\mat{A}$ we calculate
%\begin{align*}
%\mat{A} &= \frac{1}{n} \mat{X} \mat{X}^T = \frac{1}{n} \begin{pmatrix} x_{11} & \cdots & x_{1n} \\ x_{21} & \cdots & x_{2n} \end{pmatrix} \begin{pmatrix} x_{11} & \cdots & x_{1n} \\ x_{21} & \cdots & x_{2n} \end{pmatrix}^T \\
%&= \begin{pmatrix} \frac{1}{n} \sum_{i=1}^n x_{1i}^2 & \frac{1}{n} \sum_{i=1}^n x_{1i} x_{2i} \\ * & \frac{1}{n} \sum_{i=1}^n x_{2i}^2  \end{pmatrix} = \begin{pmatrix} a_{11} & a_{12} \\ * & a_{22}  \end{pmatrix}
%\end{align*}
%Note that for reasons of symmetry $a_{21} = a_{12}$, so we do not need to calculate this value. We summarize for $u,v \in \{1,2\}$
%\[ a_{uv} = \frac{1}{n} \sum_{i=1}^n x_{ui} x_{vi} \ . \]
%Analogously we get for $u,v \in \{1,2\}$
%\[ b_{uv} = \frac{1}{m} \sum_{i=1}^m y_{ui} y_{vi} \ . \]
%Now, we can calculate $\mat{C}$ according to Fig. \ref{fig:C}. Then, maximum likelihood estimation is used to obtain $\mat{M}$ and $\mat{Z}$ as described above.
%\begin{figure*}
%\begin{align*}
%\mat{C} &= \frac{1}{n \cdot m} \mat{Z} \mat{Z}^T \\
%&= \frac{1}{n \cdot m} \sum_{i=1}^n \sum_{j=1}^m \left( \begin{pmatrix} x_{1i} \\ x_{2i} \end{pmatrix} \oplus \begin{pmatrix} y_{1j} \\ y_{2j} \end{pmatrix} \right) \left( \begin{pmatrix} x_{1i} \\ x_{2i} \end{pmatrix} \oplus \begin{pmatrix} y_{1j} \\ y_{2j} \end{pmatrix} \right)^T \\
%&= \frac{1}{n \cdot m} \sum_{i=1}^n \sum_{j=1}^m \begin{pmatrix} x_{1i} y_{1j} - x_{2j} y_{2j} \\ x_{1i} y_{2j} + x_{2i} y_{1j}  \end{pmatrix} \begin{pmatrix} x_{1i} y_{1j} - x_{2j} y_{2j} \\ x_{1i} y_{2j} + x_{2i} y_{1j}  \end{pmatrix}^T \\
%&= \frac{1}{n \cdot m} \sum_{i=1}^n \sum_{j=1}^m \begin{pmatrix} x_{1i}^2 y_{1j}^2 - 2 x_{1i} x_{2i} y_{1j} y_{2j} + x_{2i}^2 y_{2j}^2  & x_{1i}^2 y_{1j} y_{2j} - x_{2i}^2 y_{1j} y_{2j} - x_{1i} x_{2i} y_{2j}^2 + x_{1i} x_{2i} y_{1j}^2 \\ * & x_{1i}^2 y_{2j}^2 + 2 x_{1i} x_{2i} y_{1j} y_{2j} + x_{2i}^2 y_{1j}^2  \end{pmatrix} \\
%&= \begin{pmatrix}
%a_{11} b_{11} - 2 a_{12} b_{12} + a_{22} b_{22} & a_{11} b_{12} - a_{22} b_{12} - a_{12} b_{22} + a_{12} b_{11} \\ * & a_{11} b_{22} + 2 a_{12} b_{12} + a_{22} b_{11}\end{pmatrix}
%\end{align*}
%\caption{Composition of two Bingham distributions. This derivation shows how to compute the scatter matrix $\mat{C}$ of the composed Bingham distribution from the scatter matrices $\mat{A}$ and $\mat{B}$ of the original Bingham distributions.}
%\label{fig:C}
%\end{figure*}
%
%Note that we can calculate $\mat{C}$ from $\mat{A}$ and $\mat{B}$ without knowing $\mat{X}$, $\mat{Y}$ and $\mat{Z}$ explicitly.  This computation can be done in an efficient way, because the solution of the equation involving the hypergeometric function is the only part which is not given in closed form. This does not present a limitation to the proposed algorithm, because there are many efficient ways for the computation of the confluent hypergeometric function of a scalar argument.

conclusion.tex
We have presented a recursive filter based on the Bingham distribution. It can be applied to circular estimation problems with $180^\circ$ symmetry. Our simulations have shown the superiority of the presented approach compared to the traditional solution of modifying a Kalman filter for the circular setting.

Future work will focus on recursive 3D pose estimation using Bingham distribution. This can be achieved by applying the presented methods in the four-dimensional case for estimating quaternions. Open challenges include an efficient estimator of the Bingham parameters based on available data. This makes an efficient evaluation of the confluent hypergeometric function necessary.


evaluation.tex
The proposed filter was evaluated in simulations. In this section, all angles are given in radians unless specified differently.

For comparison, we implemented a one-dimensional Kalman filter \cite{kalman1960}. A traditional one-dimensional Kalman filter has two issues when confronted with our situation. First, it does not take the circular nature of the problem into account. Second, it does not handle $180^\circ$ symmetry. We can circumvent both issues by restricting the estimate $x_k$ according to $0 \leq x_k \leq \pi$ and by shifting the measurement, so that $|x_k - \hat{z}_k| \leq \frac{\pi}{2}$ is satisfied.

%    %normalize Z
%    z = mod(z, pi);
%    if abs(z-Xp)>pi/2
%        z = z + pi*sign(Xp-z); 
%    end
%    
%    update
%
%    %normalize state
%    Xe = mod(Xe, pi);


In our example, we consider the estimation of an axis in robotics. This could be the axis of a symmetric rotor blade  or any other robotic joint with $180^\circ$ symmetry. We use the initial estimate with mode $(0,1)^T$
\[ \mat{M}_0^e = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}, \quad \mat{Z}_0^e = \begin{pmatrix} -1 & 0 \\ 0  & 0 \end{pmatrix} \ , \]
the system noise with mode $(1,0)^T$
\[ \mat{M}_k^w = \begin{pmatrix} 0 & 1 \\ 1 & 0 \end{pmatrix}, \quad \mat{Z}_k^w = \begin{pmatrix} -200 & 0 \\ 0  & 0 \end{pmatrix} \ , \]
and the measurement noise with mode $(1,0)^T$
\[ \mat{M}_k^v = \begin{pmatrix} 0 & 1 \\ 1 & 0 \end{pmatrix}, \quad \mat{Z}_k^v = \begin{pmatrix} -3 & 0 \\ 0  & 0 \end{pmatrix} \ . \]
The true initial state is given by $(1, 0)^T$, i.\,e., the initial estimate with mode $(0,1)^T$ is very poor. 
The initial estimate for the Kalman filter is given by 
\[ x_0^e = \atan2(\text{mode} (\mat{M}_0^e)) = \atan2(1,0) = \frac{\pi}{2} \]
and the noise means are
\[ \mu_k^w = \mu_k^v = \atan2(0,1) = 0 \ . \]

The covariance matrices for the Kalman filter are obtained by sampling the Bingham noise parameters and calculating the empirical covariance from the samples. This yields
\[ C_0^e = 0.5956, \quad C_k^w = 0.0027, \quad C_k^v=0.2836 \ , \]
which is equivalent to standard deviations of $44^\circ$ for the first time step, $3^\circ$ for the system noise and $30^\circ$ for the measurement noise.

\begin{figure*}[t]
\centering
\begin{subfigure}{0.7\textwidth}
\centering
\includegraphics[width=0.6\textwidth]{figures/2d.eps}
\caption{Ground truth and estimate.}
\end{subfigure}
\begin{subfigure}{0.7\textwidth}
\centering
\includegraphics[width=0.6\textwidth]{figures/error.eps}
\caption{Angular error.}
\end{subfigure}
\caption{An example run of Bingham and Kalman filter.}
\label{fig:examplerun}
\end{figure*}

We simulate the system for a duration of $k_{\max} = 100$ time steps. An example run is depicted in Fig. \ref{fig:examplerun}. In addition to the mode of the estimate, we plot the $95 \%$ confidence interval, which is equivalent to the $2 \sigma$ bounds in the case of the Kalman filter.

%RMSE Bingham: 0.138006
%RMSE Kalman: 0.321575

For evaluation, we consider the angular RMSE which is given by
\[ \sqrt{ \frac{1}{k_{\max}} \sum_{k=1}^{k_{\max}} (e_k)^2 } \]
with angular error 
\[ e_k = \min(\measuredangle (\vec{x}_k^\text{true}, \text{mode}(\mat{M}_k^e)) , \pi - \measuredangle (\vec{x}_k^\text{true}, \text{mode}(\mat{M}_k^e)) \]
at time step $k$. Obviously, $0 \leq e_k \leq \frac{\pi}{2}$ holds, which is consistent with our assumption of $180^\circ$ symmetry.

The presented results are based on 1000 Monte Carlo runs. Even though our filter is computationally more demanding than a Kalman filter, it is still fast enough for real-time applications. On a standard laptop, our non-optimized implementation in MATLAB needs approximately $\SI{60}{\milli\second}$ for one time step (prediction and update), which could be significantly improved by a faster evaluation of the hypergeometric function. % in Fig. \ref{fig:montecarlo-boxplot}. 
%The achievable angular RMSE in any given run strongly depends on the development of the true state and the noisy measurements, so there is a large spread in the angular RMSE. However, if we consider the ratio between the RMSE of our proposed filter and the Kalman filter, it is obvious, that our filter significantly outperforms the Kalman filter in most cases.
In Fig. \ref{fig:montecarlo-scatter}, we plot the error of our filter against the error of the Kalman filter for all runs. The proposed filter outperforms the Kalman filter in most cases, which is also true for the mean angular error in every time step as shown in Fig. \ref{fig:montecarlo-meanerror}. In particular, the significantly faster rate of convergence of the proposed filter is evident. This superiority is due to the reasons listed in the introduction. The use of the Gaussian distribution, which does not consider the problem geometry leads to suboptimal results compared to the proposed approach based on the Bingham distribution. 

In Fig. \ref{fig:montecarlo-meanerror}, we also show a comparison with a filter based on the wrapped normal distribution (denoted WN), which we previously published in \cite{kurz2013} and modified for the $180^\circ$ case. The angular error of both filters is almost indistinguishable. However, unlike the proposed filter based on the Bingham distribution, the previously published filter cannot easily be generalized to higher dimensions.

%\begin{figure*}[t]
%\centering
%\includegraphics[height=5cm]{figures/boxplot-rmse.eps}
%\includegraphics[height=5cm]{figures/boxplot-ratio.eps}
%\caption{Results of 1000 Monte Carlo runs. The angular RMSE for both filters is shown on the left, the ratio between the RMSEs in any given run is shown on the right.}
%\label{fig:montecarlo-boxplot}
%\end{figure*}

\begin{figure*}[t]
\centering
\begin{subfigure}{0.45\textwidth}
\centering
\includegraphics[height=5cm]{figures/meanerror.eps}
\caption{Time steps $1 \leq k \leq 100$.}
\end{subfigure}
\begin{subfigure}{0.45\textwidth}
\centering
\includegraphics[height=5cm]{figures/meanerror-t20.eps}
\caption{Time steps $20 \leq k \leq 100$.}
\end{subfigure}
\caption{Results of 1000 Monte Carlo runs. We show the mean error at every time step across all runs for the proposed Bingham filter, a Kalman \cite{kalman1960} filter and a filter based on the wrapped normal (WN) distribution \cite{kurz2013}. Because the initial error is large as a result of the poor initial estimate, we show two plots of different time intervals.}
\label{fig:montecarlo-meanerror}
\end{figure*}

\begin{figure}[t]
\centering
\includegraphics[width=0.6\textwidth]{figures/scatter.eps}
\caption{Results of 1000 Monte Carlo runs. Each sample represents one run. Samples below the red line indicate that the proposed filter has performed better, samples above the red line indicate that the Kalman filter has performed better.}
\label{fig:montecarlo-scatter}
\end{figure}

%    %% Parameters
%    %initial distribution
%    Z0 = [-1; 0];
%    M0 = eye(2,2);
%    %system noise
%    Zw = [-200; 0];
%    Mw =   [0 1;
%            1 0];
%	%measurement noise
%    Zv = [-3; 0];
%    Mv =   [0 1;
%            1 0];
%    %true initial state
%    x0 = [1; 0];
%    %number of time steps
%    t = 100;


filter.tex
The techniques presented in the preceding section can be applied to derive a filter based on the Bingham distribution. The system model is given by
\[ \vec{x}_{k+1} = \vec{x}_k \oplus \vec{w}_k\ ,\]
where $\vec{w}_k$ is Bingham distributed noise. The measurement model is given by
\[ \vec{z}_k = \vec{x}_k \oplus \vec{v}_k\ ,\]
where $\vec{v}_k$ is Bingham distributed noise and $\vec{x}_k$ is an uncertain Bingham distributed system state. Intuitively, this means that both system and measurement model are the identity disturbed by Bingham distributed noise. Note that $\vec{w}_k$ and $\vec{v}_k$ can include a constant offset. For example $\vec{w}_k$ could include a known angular velocity. Alternatively, to avoid dealing with biased noise distributions, a rotation may be applied to $\vec{x}_k$ first and unbiased noise added subsequently.

The predicted and estimated distributions at time $k$ are described by their parameter matrices $(\mat{M}_k^p, \mat{Z}_k^p)$ and $(\mat{M}_k^e, \mat{Z}_k^e)$ respectively. The noise distributions at time $k$ are described by $(\mat{M}_k^w, \mat{Z}_k^w)$ and $(\mat{M}_k^v,\mat{Z}_k^v)$.

\subsection{Prediction}
The prediction can be calculated according to
\[ (\mat{M}_{k+1}^p, \mat{Z}_{k+1}^p) = \text{composition}((\mat{M}_k^e, \mat{Z}_k^e),(\mat{M}_k^w, \mat{Z}_k^w)) \ , \]
which uses the previously introduced composition operation to disturb the estimate with the system noise.

\begin{algorithm}[tb]
\KwIn{estimate $\mat{M}_k^e, \mat{Z}_k^e$, noise $\mat{M}_k^w, \mat{Z}_k^w$}
\vspace{4mm}
\KwOut{prediction $\mat{M}_{k+1}^p, \mat{Z}_{k+1}^p$}
\vspace{4mm}
%composition of estimate and noise
\tcc{calculate covariance matrices $\mat{A},\mat{B}$}
$\mat{A} = \sum_{i=1}^d \frac{1}{F} \frac{\partial F}{\partial z_i} \mat{M}_k^e(:,i) \mat{M}_k^e(:,i)^T$\; %todo formular for F
$\mat{B} = \sum_{i=1}^d \frac{1}{F} \frac{\partial F}{\partial z_i} \mat{M}_k^w(:,i) \mat{M}_k^w(:,i)^T$\;
\tcc{calculate $\mat{C}$ according to Lemma \ref{lemma:composition}}
$c_{11} = a_{11} b_{11} - 2 a_{12} b_{12} + a_{22} b_{22}$\;
$c_{12} = a_{11} b_{12} - a_{22} b_{12} - a_{12} b_{22} + a_{12} b_{11}$\; 
$c_{22} = a_{11} b_{22} + 2 a_{12} b_{12} + a_{22} b_{11}$\;
$\mat{C} = \begin{pmatrix} c_{11} & c_{12} \\ c_{12} & c_{22} \end{pmatrix}$\;
\tcc{calculate $\mat{M}_{k+1}^p, \mat{Z}_{k+1}^p$ based on $\mat{C}$}
$\mat{M}_{k+1}^p, \mat{Z}_{k+1}^p \gets \text{MLE}(\mat{C})$\; %todo formula?
\label{algo:predict}
\caption{Prediction}
\end{algorithm}

\subsection{Update}
Given a measurement $\hat{\vec{z}}$, we can calculate the updated distribution according to Bayes' rule
\[ f(\mat{M}_k, \mat{Z}_k | \hat{\vec{z}}) = c \cdot f(\hat{\vec{z}} | \mat{M}_k, \mat{Z}_k) \cdot f(\mat{M}_k, \mat{Z}_k) \]
with some normalization constant $c$, which yields the update procedure
\[ (\mat{M}_k^e, \mat{Z}_k^e) = \text{multiply}((\mat{M}, \mat{Z}_k^e), (\mat{M}_k^p, \mat{Z}_k^p)) \]
with $\mat{M} = (\bar{\hat{\vec{z}}} \oplus \mat{M}_k^v) $, where $\bar{\vec{a}}$ indicates the complex conjugate of $\vec{a}$ and $\oplus$ is evaluated for each column of $\mat{M}_k^v$.


\begin{algorithm}[tb]
\KwIn{prediction $\mat{M}_k^p, \mat{Z}_k^p$, noise $\mat{M}_k^v, \mat{Z}_k^v$, measurement $\hat{\vec{z}}_k$}
\vspace{4mm}
\KwOut{estimate $\mat{M}_k^e, \mat{Z}_k^e$}
\vspace{4mm}
\tcc{rotate noise according to measurement}
$\mat{M} \gets \begin{pmatrix} \bar{\hat{\vec{z}}} \oplus \mat{M}_k^v \end{pmatrix}$\; %todo columnwise
%todo check complex conjugation
\tcc{multiply with prior distribution}
$(\mat{M}_k^e, \mat{Z}_k^e) \gets \text{multiply}((\mat{M}, \mat{Z}_k^v)), (\mat{M}_k^p, \mat{Z}_k^p))$\;
\label{algo:update}
\caption{Update}
\end{algorithm}

introduction.tex
Many estimation problems involve the task of estimating angular values. These problems include, but are not limited to, estimating the pose or orientation of objects. For example, tracking cars, ships, or airplanes may involve estimation of their current orientation or heading. Furthermore, many applications in the area of robotics or augmented reality depend on reliable estimation of the pose of certain objects. When estimating the orientation of two-way roads or relative angles of two unlabeled targets, the estimation task reduces to estimating an axis. This can be thought of as estimation of a directionless orientation or estimation with $180^\circ$ symmetry. All these estimation problems share the need for processing angular or directional data, which differs in many ways from the classical Euclidean setting. First, periodicity needs to be taken into account. This is especially important for measurement updates around $0$, respectively $2\pi$. Second, directional quantities do not lie in a vector space. Thus, there is no equivalent to a classical linear model, as there are no linear mappings. 

In many current applications, even simple estimation problems involving angular data are often considered as traditional linear or nonlinear estimation problems and handled with classical techniques such as the Kalman Filter \cite{kalman1960}, the extended Kalman Filter (EKF), or the unscented Kalman Filter (UKF) \cite{julier2004}. In a circular setting, most traditional approaches to filtering suffer from assuming a Gaussian probability density at a certain point. They fail to take into account the periodic nature of the problem and assume a linear vector space instead of a curved manifold. This shortcoming can cause poor results, in particular when the angular uncertainty is large. In certain cases, the filter may even diverge. 

Classical strategies to avoid these problems in an angular setting involve an ``intelligent'' repositioning of measurements or even discarding certain undesired measurements. Sometimes, nonlinear equality constraints have to be fulfilled, for example unit length of a vector, which makes it necessary to inflate the covariance \cite{julier2007}. There are also approaches, that use operators on a manifold to provide a local approximation of a vector space \cite{Hertzberg2013}. While these approaches yield feasible results, they still suffer from ignoring the true geometry of circular data within their probabilistic models, which are usually based on assuming a normally distributed noise. This assumption is often motivated by the Central Limit Theorem, i.e., the limit distribution of a normalized sum of i.i.d. random variables with finite variance is normally distributed \cite{shiryaev1995}. For angular data, this is not the case. Choosing a circular distribution for describing uncertainty offers possibly better results. 

In this paper, we consider the use of the Bingham distribution \cite{bingham1974} for recursive estimation of orientation. The Bingham distribution is defined on the hypersphere of arbitrary dimension, so it can be applied to problems of different dimensionality. Here, we focus on the two-dimensional case and apply our results to axis estimation. To the best of our knowledge, this is the first published attempt to create a recursive filter based on the Bingham distribution. 

The presented methods can also be applied to the four-dimensional case, which would allow the representation of unit quaternions. Unit quaternions could then be used to estimate the full 3D orientation of an object. It is well known that Quaternions avoid the singularities present in other representations such as Euler angles. Their only downsides are the fact that they must remain normalized and the property that there are two quaternions for every rotation ($q$ and $-q$). Both of these issues can elegantly be overcome by use of the Bingham distribution, since it is by definition restricted to the hypersphere and is $180^\circ$ symmetric.
%antipodally symmetric

%our contribution

This paper is structured as follows. First, we present an overview of previous work in the area of directional statistics and angular estimation (Sec. \ref{sec:relatedwork}). Then, we introduce our key idea in Sec. \ref{sec:keyidea}. In Sec. \ref{sec:binghamdistribution}, we give a detailed introduction to the Bingham distribution and we derive the necessary operations, which we will need to create a recursive Bingham filter. Based on these prerequisites, we introduce our filter in Sec. \ref{sec:filter}. We have carried out an evaluation in simulations, which is presented in Sec. \ref{sec:evaluation}. Finally, we conclude this work in Sec. \ref{sec:conclusion}.

%\begin{itemize}
%\item Many estimation problems involve the task of
%estimating angular values. This problems include,
%but are not limited to, estimating the pose or 
%orientation of objects.
%\item When angular data appears, even simple estimation
%problems are often considered as a nonlinear estimation
%problem and handled with classical techniques like the
%EKF, the UKF, (and the Particle Filter.)
%\item In a circular setting, most approaches to
%nonlinear filtering suffer from either assuming a
%Gaussian density at a certain point (i.e. in the case
%of EKF and UKF) (or being computationally demanding for
%high dimensional problems (i.e. in the case of the
%particle filter)). 
%\item issues with large uncertainties
%\item motivation for Bingham distribution: antipodally symmetric, extensible to higher dimensions, applicable to quaternions
%\item possible applications
%\item Structure of the paper.
%\end{itemize}
key_idea.tex
The goal of this paper is the derivation of a recursive filter based on the Bingham distribution. Rather than relying on the traditional Gaussian distribution, we chose to represent all occurring probability densities as Bingham. The Bingham distribution is defined on the hypersphere and is antipodally symmetric, which makes it interesting for applications in angular estimation with inherent $180^\circ$ symmetry and for problems where $180^\circ$ symmetry occurs as a result of parameterization, e.g., in the case of quaternions. Although we restrict ourselves to the two-dimensional case in this paper, we would like to emphasize that most of presented methods are easily generalized to higher dimensions.

In order to derive a recursive filter, we need to be able to perform two operations. First, we need to calculate the predicted state at the next time step from the current state and the system noise affecting the state. In a traditional estimation problem in $\mathbb{R}^d$ with additive noise, this involves a convolution with the noise density. We provide a suitable analogue on the hypersphere, which we call \emph{composition}. Since Bingham distributions are not closed under compositions, we present an approximate solution to this problem, which is based on matching covariance matrices. 

Second, we need to perform a Bayes update. As usual, this requires the \emph{multiplication} of the prior density with the likelihood density. We prove that Bingham distributions are closed under multiplication and show how to obtain the posterior density.



main.tex
\documentclass[preprint,1p,11pt]{IR-Template/ISAS_IR}

\usepackage[cmex10]{amsmath}
\usepackage{amsfonts}
\usepackage[ruled]{algorithm2e}
\usepackage[english]{babel}
%\usepackage[utf8]{inputenc}
%\usepackage{packages-ig/abkuerzungen}
\usepackage{ISASPackages/ISASMatheMacros}
\usepackage{graphicx}
\usepackage{epstopdf}
\usepackage{subcaption}
\usepackage{flushend}
\usepackage{siunitx}
\usepackage{tabularx}
\usepackage{hyperref}

\usepackage{amsthm}
\newtheorem{definition}{Definition}
\newtheorem{lemma}{Lemma}


\def\PDF#1{#1}

%\input{IR-Template/Defs.h} % Lemma, Theorem etc.

\begin{document}

\begin{frontmatter}

\title{Recursive Estimation of Orientation \\ Based on the Bingham Distribution}
\tnotetext[t1]{Draft submitted to 16th International Conference on Information FUSION on March 15, 2013}

\author[isas]{Gerhard~Kurz}
\ead{gerhard.kurz@kit.edu}

\author[isas]{Igor~Gilitschenski}
\ead{gilitschenski@kit.edu}

\author[ucl]{Simon~Julier}
\ead{s.julier@cs.ucl.ac.uk}

\author[isas]{Uwe~D.~Hanebeck}
\ead{uwe.hanebeck@ieee.org}

\address[isas]{Intelligent Sensor-Actuator-Systems Laboratory (ISAS)\\
Institute for Anthropomatics\\
Karlsruhe Institute of Technology (KIT), Germany\vspace{3mm}}

\address[ucl]{Virtual Environments and Computer Graphics Group\\
Department of Computer Science\\
University College London (UCL), United Kingdom}

\begin{abstract}
Directional estimation is a common problem in many tracking applications. Traditional filters such as the Kalman filter perform poorly because they fail to take the periodic nature of the problem into account. We present a recursive filter for directional data based on the Bingham distribution in two dimensions. The proposed filter can be applied to circular filtering problems with 180 degree symmetry, i.e., rotations by 180 degrees cannot be distinguished. It is easily implemented using standard numerical techniques and suitable for real-time applications. The presented approach is extensible to quaternions, which allow tracking arbitrary three-dimensional orientations. We evaluate our filter in a challenging scenario and compare it to a traditional Kalman filtering approach.
\end{abstract}

\end{frontmatter}

\section{Introduction} \label{sec:introduction}
\noindent
\input{introduction}

\section{Related Work} \label{sec:relatedwork}
\noindent
\input{related_work}

\section{Key Idea of the Bingham Filter} \label{sec:keyidea}
\noindent
\input{key_idea}

\section{Bingham Distribution} \label{sec:binghamdistribution}
\noindent
\input{bingham_distribution}

\section{Filter Implementation} \label{sec:filter}
\noindent
\input{filter}

\section{Evaluation} \label{sec:evaluation}
\noindent
\input{evaluation}

\section{Conclusion} \label{sec:conclusion}
\noindent
\input{conclusion}

\section*{Acknowledgment} \noindent
This work was partially supported by grants from the German Research Foundation
(DFG) within the Research Training Groups RTG 1194 ``Self-organizing
Sensor-Actuator-Networks'' and RTG 1126 ``Soft-tissue Surgery: New
Computer-based Methods for the Future Workplace''.


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    
\begin{appendix}


\input{appendix.tex}
	
\end{appendix}
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\bibliographystyle{ieeetemplates/IEEEtran}
\bibliography{literature-ig,gk-bingham}


%\printbibliography 
%\bibliographystyle{IEEEtran}
%\bibliography{%
 	%ISASPublikationen.bib%
%}

\end{document}
related_work.tex
Directional statistics is a subdiscipline of statistics, which focuses on dealing with directional data. Classical results in directional statistics are summed up in the books by Mardia and Jupp \cite{mardia1999} and by Jammalamadaka and Sengupta \cite{jammalamadaka2001}. Directional statistics differs from traditional statistics by the fact that random variables located on manifolds (for example the circle or the sphere) are considered rather than random variables located in vector spaces (typically $\mathbb{R}^d$).

There is a broad range of research for investigating the 2D orientation, e.\,g., the work by Krindis et al. \cite{krinidis2006}. A recursive filter based on the von Mises distribution for estimating the orientation on the $SO(2)$ was presented in \cite{azmani2009}. Later, a nonlinear filter based on von Mises and wrapped normal distributions was presented in \cite{kurz2013}.

In 1974, Bingham proposed his distribution in \cite{bingham1974}. Further work on the Bingham distribution has been done by Kent \cite{kent1987} as well as Jupp and Mardia \cite{jupp1979}. So far, there have only been a few applications of the Bingham distribution, for example in geology \cite{kunze2004}. In 2011, Glover used the Bingham distribution for a Monte Carlo based pose estimation \cite{glover2011}. Glover also released a library called libbingham \cite{Glover13} that includes implementations of some of the methods discussed in Sec. \ref{sec:binghamdistribution}. It should be noted that our implementation is not based on libbingham.

%hypergeom stuff? quaternions?

%\begin{itemize}
%  \item Classical results in directional statistics are summed up in 
%  \cite{mardia1999,jammalamadaka2001}. 
%  \item Bingham proposed his distribution in \cite{bingham1974}. 
%  \item A nonlinear filter for estimating the orientation on the SO(2) was presented 
%  in \cite{kurz2013}. Furthermore another recursive filter was presented in 
%  \cite{azmani2009}.
%  \item In \cite{glover2011}, Glover used the Bingham distribution for a Monte Carlo based pose
%  estimation.
%  \item A broad range of research for investigating the 2D orientation exists, 
%  e.g.  \cite{krinidis2006}.
%\end{itemize}
main.bbl
% Generated by IEEEtran.bst, version: 1.12 (2007/01/11)
\begin{thebibliography}{10}
\providecommand{\url}[1]{#1}
\csname url@samestyle\endcsname
\providecommand{\newblock}{\relax}
\providecommand{\bibinfo}[2]{#2}
\providecommand{\BIBentrySTDinterwordspacing}{\spaceskip=0pt\relax}
\providecommand{\BIBentryALTinterwordstretchfactor}{4}
\providecommand{\BIBentryALTinterwordspacing}{\spaceskip=\fontdimen2\font plus
\BIBentryALTinterwordstretchfactor\fontdimen3\font minus
  \fontdimen4\font\relax}
\providecommand{\BIBforeignlanguage}[2]{{%
\expandafter\ifx\csname l@#1\endcsname\relax
\typeout{** WARNING: IEEEtran.bst: No hyphenation pattern has been}%
\typeout{** loaded for the language `#1'. Using the pattern for}%
\typeout{** the default language instead.}%
\else
\language=\csname l@#1\endcsname
\fi
#2}}
\providecommand{\BIBdecl}{\relax}
\BIBdecl

\bibitem{kalman1960}
R.~E. Kalman, ``{A} {N}ew {A}pproach to {L}inear {F}iltering and {P}rediction
  {P}roblems,'' \emph{Transactions of the {ASME} Journal of Basic Engineering},
  vol.~82, pp. 35--45, 1960.

\bibitem{julier2004}
S.~Julier and J.~Uhlmann, ``Unscented filtering and nonlinear estimation,''
  \emph{Proceedings of the IEEE}, vol.~92, no.~3, pp. 401--422, mar 2004.

\bibitem{julier2007}
S.~Julier and J.~LaViola, ``On kalman filtering with nonlinear equality
  constraints,'' \emph{IEEE Transactions on Signal Processing}, vol.~55, no.~6,
  pp. 2774--2784, 2007.

\bibitem{Hertzberg2013}
C.~Hertzberg, R.~Wagner, U.~Frese, and L.~Schr�der, ``Integrating generic
  sensor fusion algorithms with sound state representations through
  encapsulation of manifolds,'' \emph{Information Fusion}, vol.~14, no.~1, pp.
  57 -- 77, 2013.

\bibitem{shiryaev1995}
A.~N. Shiryaev, \emph{Probability}, 2nd~ed.\hskip 1em plus 0.5em minus
  0.4em\relax Springer, 1995.

\bibitem{bingham1974}
C.~Bingham, ``An antipodally symmetric distribution on the sphere,'' \emph{The
  Annals of Statistics}, vol.~2, no.~6, pp. 1201--1225, Nov. 1974.

\bibitem{mardia1999}
K.~V. Mardia and P.~E. Jupp, \emph{Directional Statistics}, 1st~ed.\hskip 1em
  plus 0.5em minus 0.4em\relax Wiley, 1999.

\bibitem{jammalamadaka2001}
S.~R. Jammalamadaka and A.~Sengupta, \emph{Topics in Circular
  Statistics}.\hskip 1em plus 0.5em minus 0.4em\relax World Scientific Pub Co
  Inc, 2001.

\bibitem{krinidis2006}
S.~Krinidis and V.~Chatzis, ``Frequency-based object orientation and scaling
  determination,'' in \emph{2006 {IEEE} International Symposium on Circuits and
  Systems, 2006. {ISCAS} 2006. Proceedings}, May 2006, p. 4 pp.

\bibitem{azmani2009}
M.~Azmani, S.~Reboul, J.-B. Choquel, and M.~Benjelloun, ``A recursive fusion
  filter for angular data,'' in \emph{2009 {IEEE} International Conference on
  Robotics and Biomimetics ({ROBIO)}}, Dec. 2009, pp. 882 --887.

\bibitem{kurz2013}
G.~Kurz, I.~Gilitschenski, and U.~D. Hanebeck, ``Recursive nonlinear filtering
  for angular data based on circular distributions,'' in \emph{Proceedings of
  the 2013 American Control Conference ({ACC} 2013) (to appear)}, Washington
  {D.C.}, {USA}, Jun. 2013.

\bibitem{kent1987}
J.~T. Kent, ``Asymptotic expansions for the bingham distribution,''
  \emph{Journal of the Royal Statistical Society. Series C (Applied
  Statistics)}, vol. 36 (2), pp. 139--144, 1987.

\bibitem{jupp1979}
P.~E. Jupp and K.~V. Mardia, ``Maximum likelihood estimators for the matrix von
  mises-fisher and bingham distributions,'' \emph{Annals of Statistics}, vol. 7
  (3), pp. 599--606, 1979.

\bibitem{kunze2004}
K.~Kunze and H.~Schaeben, ``\BIBforeignlanguage{English}{The bingham
  distribution of quaternions and its spherical radon transform in texture
  analysis},'' \emph{\BIBforeignlanguage{English}{Mathematical Geology}},
  vol.~36, pp. 917--943, 2004.

\bibitem{glover2011}
J.~Glover, R.~Rusu, and G.~Bradski, ``Monte carlo pose estimation with
  quaternion kernels and the bingham distribution,'' in \emph{Proceedings of
  Robotics: Science and Systems}, Los Angeles, {CA}, {USA}, Jun. 2011.

\bibitem{Glover13}
\BIBentryALTinterwordspacing
J.~Glover, ``libbingham bingham statistics library,'' 2013. [Online].
  Available: \url{http://code.google.com/p/bingham/}
\BIBentrySTDinterwordspacing

\bibitem{kume2005}
A.~Kume and A.~T.~A. Wood, ``Saddlepoint approximations for the bingham and
  {Fisher–Bingham} normalising constants,'' \emph{Biometrika}, vol.~92,
  no.~2, pp. 465--476, 2005.

\bibitem{herz1955}
C.~S. Herz, ``Bessel functions of matrix argument,'' \emph{Annals of
  Mathematics}, vol.~61, no.~3, pp. 474--523, 1955.

\bibitem{koev2006}
P.~Koev and A.~Edelman, ``The efficient evaluation of the hypergeometric
  function of a matrix argument,'' \emph{Math. Comp.}, vol.~75, pp. 833--846,
  2006.

\bibitem{muirhead1982}
R.~Muirhead, \emph{Aspects of multivariate statistical theory}, ser. Wiley
  series in probability and mathematical statistics: Probability and
  mathematical statistics.\hskip 1em plus 0.5em minus 0.4em\relax Wiley, 1982.

\bibitem{abramowitz1964}
M.~{Abramowitz} and I.~A. {Stegun}, \emph{Handbook of Mathematical Functions
  with Formulas, Graphs, and Mathematical Tables}, ninth dover printing, tenth
  gpo printing~ed.\hskip 1em plus 0.5em minus 0.4em\relax New York: Dover,
  1964.

\bibitem{luke1977}
Y.~Luke, \emph{Algorithms for the computation of mathematical functions}, ser.
  Computer science and applied mathematics.\hskip 1em plus 0.5em minus
  0.4em\relax Academic Press, 1977.

\bibitem{muller2001}
K.~E. Muller, ``\BIBforeignlanguage{en}{Computing the confluent hypergeometric
  function, m(a,b,x)},'' \emph{\BIBforeignlanguage{en}{Numerische Mathematik}},
  vol.~90, no.~1, pp. 179--196, 2001.

\end{thebibliography}
