function p = getPdf_VM(alpha, mean, kappa)
%GETPDF_VM Probability in von Mises distribution.
%   P = GETPDF_VM(ALPHA, MEAN, KAPPA) returns the probability density of
%   the angle ALPHA (in radians) in a von Mises distribution of preferred
%   direction MEAN (default is 0) and concentration kappa (default is 1).
%
%   The parameter ALPHA can also contain several angles (dimension 1xN),
%   P is then of dimension 1xN.

%   Author: Damien Teney

% Check the arguments
if nargin < 3
  kappa = 1;
end
if nargin < 2
  mean = 0;
end

% Make sure we have 1xN vectors
alpha = alpha(:)';
mean = mean(:)';

% Make sure kappa is an integer
kappa = round(kappa);

% Get the normalization constant
%normalizationConstant = 1 / (2 * pi * besseli(0, kappa));
% Precomputed values (1:512)
normalizationConstants = [1.257083e-001, 6.981750e-002, 3.260842e-002, 1.408211e-002, 5.842720e-003, 2.367165e-003, 9.440136e-004, 3.722364e-004, 1.455346e-004, 5.652378e-005, 2.183648e-005, 8.399154e-006, 3.218861e-006, 1.229769e-006, 4.685860e-007, 1.781360e-007, 6.758257e-008, 2.559414e-008, 9.677314e-009, 3.653839e-009, 1.377798e-009, 5.189390e-010, 1.952481e-010, 7.338999e-011, 2.756139e-011, 1.034214e-011, 3.877855e-012, 1.453011e-012, 5.440811e-013, 2.036083e-013, 7.615204e-014, 2.846674e-014, 1.063600e-014, 3.972064e-015, 1.482734e-015, 5.532615e-016, 2.063609e-016, 7.694203e-017, 2.867789e-017, 1.068529e-017, 3.980041e-018, 1.482034e-018, 5.517013e-019, 2.053199e-019, 7.639141e-020, 2.841512e-020, 1.056698e-020, 3.928733e-021, 1.460357e-021, 5.427179e-022, 2.016515e-022, 7.491080e-023, 2.782315e-023, 1.033214e-023, 3.836178e-024, 1.424082e-024, 5.285683e-025, 1.961552e-025, 7.278361e-026, 2.700252e-026, 1.001646e-026, 3.715055e-027, 1.377715e-027, 5.108557e-028, 1.894016e-028, 7.021297e-029, 2.602560e-029, 9.645738e-030, 3.574562e-030, 1.324537e-030, 4.907507e-031, 1.818086e-031, 6.734812e-032, 2.494570e-032, 9.239018e-033, 3.421505e-033, 1.266983e-033, 4.691236e-034, 1.736873e-034, 6.430039e-035, 2.380264e-035, 8.810557e-036, 3.260987e-036, 1.206877e-036, 4.466282e-037, 1.652719e-037, 6.115362e-038, 2.262646e-038, 8.371104e-039, 3.096859e-039, 1.145600e-039, 4.237584e-040, 1.567393e-040, 5.797117e-041, 2.143984e-041, 7.928791e-042, 2.932031e-042, 1.084194e-042, 4.008878e-043, 1.482232e-043, 5.480093e-044, 2.025994e-044, 7.489751e-045, 2.768701e-045, 1.023445e-045, 3.782974e-046, 1.398243e-046, 5.167886e-047, 1.909961e-047, 7.058584e-048, 2.608512e-048, 9.639405e-049, 3.561970e-049, 1.316174e-049, 4.863170e-050, 1.796839e-050, 6.638694e-051, 2.452676e-051, 9.061124e-052, 3.347407e-052, 1.236573e-052, 4.567898e-053, 1.687323e-053, 6.232547e-054, 2.302071e-054, 8.502723e-055, 3.140390e-055, 1.159833e-055, 4.283456e-056, 1.581903e-056, 5.841879e-057, 2.157310e-057, 7.966363e-058, 2.941679e-058, 1.086221e-058, 4.010785e-059, 1.480910e-059, 5.467846e-060, 2.018796e-060, 7.453452e-061, 2.751765e-061, 1.015907e-061, 3.750474e-062, 1.384547e-062, 5.111148e-063, 1.886770e-063, 6.964811e-064, 2.570926e-064, 9.489860e-065, 3.502840e-065, 1.292918e-065, 4.772130e-066, 1.761344e-066, 6.500796e-067, 2.399274e-067, 8.854908e-068, 3.267980e-068, 1.206051e-068, 4.450856e-069, 1.642527e-069, 6.061403e-070, 2.236791e-070, 8.254091e-071, 3.045825e-071, 1.123913e-071, 4.147172e-072, 1.530255e-072, 5.646348e-073, 2.083357e-073, 7.686919e-074, 2.836177e-074, 1.046422e-074, 3.860762e-075, 1.424400e-075, 5.255132e-076, 1.938779e-076, 7.152631e-077, 2.638739e-077, 9.734647e-078, 3.591180e-078, 1.324791e-078, 4.887097e-079, 1.802802e-079, 6.650258e-080, 2.453141e-080, 9.048993e-081, 3.337887e-081, 1.231223e-081, 4.541464e-082, 1.675131e-082, 6.178680e-083, 2.278960e-083, 8.405658e-084, 3.100280e-084, 1.143469e-084, 4.217373e-085, 1.555443e-085, 5.736677e-086, 2.115735e-086, 7.802909e-087, 2.877706e-087, 1.061282e-087, 3.913904e-088, 1.443391e-088, 5.322954e-089, 1.962982e-089, 7.238934e-090, 2.669487e-090, 9.844102e-091, 3.630107e-091, 1.338621e-091, 4.936182e-092, 1.820203e-092, 6.711870e-093, 2.474928e-093, 9.125926e-094, 3.365012e-094, 1.240771e-094, 4.575011e-095, 1.686895e-095, 6.219844e-096, 2.293330e-096, 8.455690e-097, 3.117649e-097, 1.149479e-097, 4.238094e-098, 1.562557e-098, 5.760988e-099, 2.123997e-099, 7.830808e-100, 2.887056e-100, 1.064387e-100, 3.924102e-101, 1.446694e-101, 5.333463e-102, 1.966246e-102, 7.248739e-103, 2.672288e-103, 9.851451e-104, 3.631728e-104, 1.338822e-104, 4.935468e-105, 1.819409e-105, 6.707002e-106, 2.472424e-106, 9.114103e-107, 3.359706e-107, 1.238469e-107, 4.565258e-108, 1.682837e-108, 6.203192e-109, 2.286573e-109, 8.428522e-110, 3.106808e-110, 1.145181e-110, 4.221144e-111, 1.555905e-111, 5.734991e-112, 2.113874e-112, 7.791521e-113, 2.871853e-113, 1.058519e-113, 3.901507e-114, 1.438013e-114, 5.300173e-115, 1.953504e-115, 7.200047e-116, 2.653709e-116, 9.780665e-117, 3.604794e-117, 1.328585e-117, 4.896610e-118, 1.804673e-118, 6.651182e-119, 2.451298e-119, 9.034218e-120, 3.329524e-120, 1.227074e-120, 4.522273e-121, 1.666633e-121, 6.142148e-122, 2.263591e-122, 8.342049e-123, 3.074290e-123, 1.132959e-123, 4.175237e-124, 1.538669e-124, 5.670309e-125, 2.089611e-125, 7.700551e-126, 2.837759e-126, 1.045747e-126, 3.853675e-127, 1.420107e-127, 5.233167e-128, 1.928437e-128, 7.106309e-129, 2.618666e-129, 9.649698e-130, 3.555862e-130, 1.310309e-130, 4.828363e-131, 1.779196e-131, 6.556099e-132, 2.415821e-132, 8.901878e-133, 3.280169e-133, 1.208672e-133, 4.453674e-134, 1.641066e-134, 6.046883e-135, 2.228100e-135, 8.209859e-136, 3.025063e-136, 1.114631e-136, 4.107005e-137, 1.513273e-137, 5.575803e-138, 2.054448e-138, 7.569739e-139, 2.789103e-139, 1.027652e-139, 3.786390e-140, 1.395091e-140, 5.140174e-141, 1.893874e-141, 6.977857e-142, 2.570935e-142, 9.472360e-143, 3.489983e-143, 1.285838e-143, 4.737480e-144, 1.745446e-144, 6.430779e-145, 2.369293e-145, 8.729150e-146, 3.216054e-146, 1.184875e-146, 4.365361e-147, 1.608295e-147, 5.925288e-148, 2.182988e-148, 8.042505e-149, 2.962985e-149, 1.091606e-149, 4.021613e-150, 1.481607e-150, 5.458380e-151, 2.010911e-151, 7.408326e-152, 2.729264e-152, 1.005470e-152, 3.704173e-153, 1.364619e-153, 5.027245e-154, 1.852025e-154, 6.822788e-155, 2.513479e-155, 9.259486e-156, 3.411119e-156, 1.256623e-156, 4.629262e-157, 1.705362e-157, 6.282319e-158, 2.314311e-158, 8.525536e-159, 3.140654e-159, 1.156956e-159, 4.261987e-160, 1.570022e-160, 5.783593e-161, 2.130532e-161, 7.848325e-162, 2.891108e-162, 1.065001e-162, 3.923144e-163, 1.445163e-163, 5.323509e-164, 1.961000e-164, 7.223633e-165, 2.660923e-165, 9.801834e-166, 3.610613e-166, 1.330004e-166, 4.899184e-167, 1.804650e-167, 6.647534e-168, 2.448651e-168, 9.019691e-169, 3.322424e-169, 1.223818e-169, 4.507934e-170, 1.660491e-170, 6.116378e-171, 2.252945e-171, 8.298614e-172, 3.056744e-172, 1.125930e-172, 4.147267e-173, 1.527606e-173, 5.626774e-174, 2.072555e-174, 7.633989e-175, 2.811872e-175, 1.035710e-175, 3.814868e-176, 1.405140e-176, 5.175571e-177, 1.906320e-177, 7.021531e-178, 2.586227e-178, 9.525773e-179, 3.508589e-179, 1.292300e-179, 4.759849e-180, 1.753160e-180, 6.457268e-181, 2.378344e-181, 8.759906e-182, 3.226435e-182, 1.188352e-182, 4.376895e-183, 1.612077e-183, 5.937512e-184, 2.186864e-184, 8.054489e-185, 2.966559e-185, 1.092614e-185, 4.024197e-186, 1.482145e-186, 5.458847e-187, 2.010527e-187, 7.404877e-188, 2.727248e-188, 1.004454e-188, 3.699428e-189, 1.362504e-189, 5.018108e-190, 1.848166e-190, 6.806766e-191, 2.506915e-191, 9.232881e-192, 3.400429e-192, 1.252360e-192, 4.612363e-193, 1.698700e-193, 6.256174e-194, 2.304092e-194, 8.485742e-195, 3.125206e-195, 1.150977e-195, 4.238901e-196, 1.561129e-196, 5.749413e-197, 2.117420e-197, 7.798114e-198, 2.871911e-198, 1.057673e-198, 3.895209e-199, 1.434528e-199, 5.283071e-200, 1.945641e-200, 7.165361e-201, 2.638836e-201, 9.718194e-202, 3.578968e-202, 1.318042e-202, 4.853995e-203, 1.787593e-203, 6.583197e-204, 2.424399e-204, 8.928335e-205, 3.288030e-205, 1.210877e-205, 4.459268e-206, 1.642200e-206, 6.047661e-207, 2.227142e-207, 8.201769e-208, 3.020412e-208, 1.112305e-208, 4.096194e-209, 1.508469e-209, 5.555092e-210, 2.045715e-210, 7.533524e-211, 2.774280e-211, 1.021648e-211, 3.762285e-212, 1.385483e-212, 5.102107e-213, 1.878871e-213, 6.919006e-214, 2.547942e-214, 9.382841e-215, 3.455241e-215, 1.272393e-215, 4.685583e-216, 1.725460e-216, 6.353974e-217, 2.339834e-217, 8.616358e-218, 3.172938e-218, 1.168419e-218, 4.302638e-219, 1.584419e-219, 5.834512e-220, 2.148514e-220, 7.911720e-221, 2.913418e-221, 1.072837e-221, 3.950610e-222];
normalizationConstant = normalizationConstants(kappa);

% Evaluate pdf
p = normalizationConstant * exp(kappa * cos(alpha - mean));
