Erstellen einer Spiralstruktur in Python, mit Hyperbolic Tangens

Question

Ich versuche, eine spiralförmige Struktur, wie Spiralarme einer Galaxie, in einem 2D-Array in Python zu erstellen. Die erste und einfache Art, wie ich es tat, war eine einfache Log-Spiral-Funktion, wie im Bild definiert: log spiral function

Die x und y Werte von

x,y=meshgrid(arange(0,M=400,1), arange(0,N=400,1)) 

M und N erstellt werden, sind die Dimensionen des Arrays. Der Radius Koordinate ist einfach, wie die Gleichung des letzten Bildes,

r=(abs(x-gal_center[1])**(2.0)+((abs(y-gal_center[0]))/(q))**(2.0))**(0.5) 

das Profil Helligkeit von f (r) zu schaffen und ploting

plt.imshow((abs(galaxy_model))**0.2) 

gib mir eine commom spiralförmigen Struktur, wie eine Spiral Galaxis.

Eine andere Möglichkeit, dies zu tun, ist die Verwendung einer anderen Funktion, die hyperbolic tangent. In den Gleichungen des letzten Bildes, außer r, das wie zuvor definiert ist, sind alle anderen Parameter einstellbare Zahlen.

Für diese Funktion habe ich Probleme, eine Spiralstruktur in einem 2D-Array zu machen. Ich weiß nicht, ob ich den hyperbolischen Tangens verwenden muss, um eine Koordinatentransformation im Array oder eine Matrix-/Array-Verzerrung vorzunehmen, um eine Spiralstruktur zu erzeugen. Ich habe es versucht, aber ich konnte es nicht.

Wie kann ich verfahren, um dieses spira/image zu erstellen, mit den oben genannten Definitionen? Danke für die Hilfe!

Weitere Informationen zu diesem Thema, in den Referenzen:

1. Peng, Y. Chien et al; Detaillierte strukturelle Zerlegung von Galaxienbildern, 2002
2. Peng, Y. Chien et al .; Detaillierte Zerlegung von Galaxy-Bildern. II. Jenseits achsensymmetrischer Modelle, 2009
3. Peng, Y. Chien, Galfit Benutzerhandbuch, 2003
4. Rowe, Barnaby et al; GALSIM: Die modulare Galaxie Bildsimulation Toolkit 2015

Edited:

Der Code, den ich verwende, wie folgt:

from __future__ import division 
import numpy as np 
from numpy import* 
import matplotlib.pyplot as pyplot 
import scipy as sp 
from scipy import* 
import pylab as pl 
from pylab import* 
import math 
from math import* 
import pyfits as pf 
from pyfits import* 

def exponential_profile(Io,ro,r): 
    Iexp=0.5*Io*np.exp(-r/ro) 
    return Iexp 

def sersic_profile(Io,ro,r,n): 
    Iser=Io*np.exp(-(r/ro)**(1/n)) 
    return Iser 

def galaxy_model1(q,c,gal_center,Io,ro,n,M,N,xi,p,n1,n2,s1,s2,k): 
    x,y=meshgrid(arange(-M/2,M/2,1), arange(-N/2,N/2,1)) 
    r=(abs(x-0*gal_center[1])**(c+2.0)+((abs(y-0*gal_center[0]))/(q))**(c+2.0))**(1.0/(c+2.0)) 
    power=2.0 
    fr=(30-xi*np.log(1.0+r**power)+(1.0/p)*np.cos(n1*arctan2(x,y)+k*np.log(s1+r**power))+(1.0/p)*np.cos(n2*arctan2(x,y)+k*np.log(s2+r**power)) ) 
    I_exp=exponential_profile(Io,ro,r) 
    I_ser=sersic_profile(Io,ro,r,n) 
    galaxy_model_1=0.1*I_exp+0.1*I_ser+0.5*fr 
    return galaxy_model_1 

def galaxy_model2(q,c,Cb,rout,rin,Oout,a,M,N,Io,ro,n): 
    gal_center=(M/2,N/2) 
    x,y=meshgrid(arange(0,M,1), arange(0,N,1)) 
    r=(abs(x-0*gal_center[1])**(c+2.0)+((abs(y-0*gal_center[0]))/(q))**(c+2.0))**(1.0/(c+2.0)) 
    A=2*Cb/(abs(Oout)+Cb)-1.00001 
    B=(2-np.arctanh(A))*((rout)/(rout-rin)) 
    T=0.5*(np.tanh(B*(r/rout-1)+2)+1) 
    Or=Oout*T*(0.5*(r/rout+1))**a 
    I_exp=exponential_profile(Io,ro,r) 
    I_ser=sersic_profile(Io,ro,r,n) 
    galaxy_model_2=0.1*I_exp+0.1*I_ser+0.5*Or 
    return galaxy_model_2 
galaxy_model_1=galaxy_model1(q,c,(M/2,N/2),Io,ro,n,M,N,xi,p,n1,n2,s1,s2,k) 
galaxy_model_2=galaxy_model2(q,c,Cb,rout,rin,Oout,a,M,N,Io,ro,n) 
fig=plt.figure() 
ax1=fig.add_subplot(121) 
ax1.imshow((abs(galaxy_model_1))**0.2) 
pf.writeto('gal_1.fits', galaxy_model_1, clobber=1) 
ax2=fig.add_subplot(122, axisbg='white') 
ax2.imshow((abs(galaxy_model_2))**0.2) 
plt.show() 

Eine Reihe von Parametern sein:

M=400 
N=400 
q=0.8 
c=0.0 
Io=100.0 
ro=10.0 
n=3.0 
xi=2.0 
p=1.7 
n1=3.0 
n2=3.0 
s1=0.05 
s2=0.5 
k=3.0 
Cb=0.23 
rout=100.0 
rin=10.0 
Oout=pi/2 
a=0.0 

Answer 1

Ich bin mir nicht sicher, dass das genau richtig ist, aber ich denke, es ist nahe, und produziert ähnliche Ergebnisse wie das Papier:

def galaxy_model2(q,c,Cb,rout,rin,Oout,a,M,N,Io,ro,n): 
    gal_center=(0,0) 
    x,y=meshgrid(arange(-M/2,M/2,1), arange(-N/2,N/2,1)) 
    r=(abs(x-gal_center[1])**(c+2.0)+((abs(y-gal_center[0]))/(q))**(c+2.0))**(1.0/(c+2.0)) 
    A=2*Cb/(abs(Oout)+Cb)-1.00001 
    B=(2-np.arctanh(A))*((rout)/(rout-rin)) 
    T=0.5*(np.tanh(B*(r/rout-1)+2)+1) 
    Or=Oout*T*(0.5*(r/rout+1))**a 
    Or=30-np.log(1.0+r**2.0)+(2.0/p)*np.cos(n2*arctan2(x,y)+k*Or) 
    I_exp=exponential_profile(Io,ro,r) 
    I_ser=sersic_profile(Io,ro,r,n) 
    #galaxy_model_2=0.5*Or 
    return Or 

Die einzige Änderung ist, dass ich

Or=30-np.log(1.0+r**2.0)+(2.0/p)*np.cos(n2*arctan2(x,y)+k*Or) 

verwenden, um eine Galaxie Plot zu erstellen.

np.cos(n1*arctan2(x,y)) 

schafft dieses Grundstück:

und ich kann es drehen sich durch k*Or

diese Verwendung mit n2 = 3 Hinzufügen erhalte ich: