Мне нужно локализовать событие на небе и создать карту с помощью heelpy.
Во время тестовой сессии я использовал метод фон Мизеса-Фишера, так как предположил, что дисперсия для $'\theta'$ такая же, как и для $'\phi'$. Все работало хорошо, и мне удалось получить параметр концентрации $'k'$, используя $'k=1/\sigma^2'$.
Моя функция для оценки вероятности в пиксель был
Код: Выделить всё
def Mises_Fisher(theta,phi,DS_theta,DS_phi,conc):
meanvec=hp.ang2vec(DS_theta,DS_phi)
meanvec=np.asarray(meanvec,dtype=np.float128)
norm=np.sqrt(np.dot(meanvec,meanvec))
meanvec=meanvec/norm
var=hp.ang2vec(theta,phi)
var=np.asarray(var,dtype=np.float128)
norm=np.sqrt(np.dot(var,var))
var=var/norm
factor=np.dot(conc*var,meanvec)
factor=np.float128(factor)
#Normalization is futile, we will devide by the sum
#fullnorm=conc/(2*np.pi*(np.exp(conc)-np.exp(-conc)))
ret=np.float128(np.exp(factor))#/fullnorm
#ret=factor
return ret
Теперь у меня есть ковариационная матрица, которая больше не является изотропной, как я могу генерировать точки на сфере в этом случае? Я знаю о существовании распределения Кента, но у меня нет необходимых величин, таких как $'\gamma'$,$'\beta'$ и $'k'$. Я не знаю, есть ли способ получить их с помощью ковариационной матрицы.
Чтобы еще больше усложнить проблему, матрица имеет 11 измерений, но мне нужна только проекция на ( расстояние, тета, фи)-пространство, чтобы затем создать карту.
На данный момент я попытался извлечь из многомерного гауссиана, но поскольку гауссиана находится в касательной плоскости, для некоторых отклонений распределение создают углы, выходящие за пределы области определения для $'\theta'$ и $'\phi'$.
Вот код, который я использую для извлечения точек
Код: Выделить всё
num_samples = 10**7
samples = np.random.multivariate_normal(perm_mean, perm_cov, num_samples)
phi = samples[:, 2]
theta = samples[:, 1]
print(np.min(theta),np.max(theta))
print(np.min(phi),np.max(phi))
print('starting mean values')
print('theta={}, phi={}'.format(perm_mean[1],perm_mean[2]))
Заранее всем спасибо!
Подробнее здесь: https://stackoverflow.com/questions/790 ... ntration-p