Pixel-zu-RA / DEC-Kartierung in der digitalisierten Astrofotografie

Ich habe ein 1443x998-Bild der Sterne (aufgenommen mit einer 35-mm-Kamera und dann gescannt) mit den folgenden Sternen an den folgenden Pixelpositionen:

Altair x=782, y=532 [19h50m46.9990s RA, +08 52'05.959'' DEC] 
Sualocin, x=311, y=146 [20h 39m 38.287s +15 54'43.49'' DEC] 
Denebokab, x=1023, y=815 [19h25m29.9005s +03 06' 53.191'' DEC] 

Welche mathematische Funktion konvertiert die Pixelposition in RA / DEC und umgekehrt? Anmerkungen:

• Helle Sterne sind Blobs auf dem Bild. Die obigen Koordinaten sind ungefähr die Mitte des Blobs, können jedoch um + -2 Pixel abweichen.

• Ich weiß, dass ich die Himmelskugel so drehen kann, dass die Bildmitte Polarkoordinaten 0,0 hat. Die eigentliche Frage ist also "wie man diese Rotation findet" (aber siehe nächster Punkt).

• Wenn Höhe / Azimut in Bildern linear wären, wäre dies einfach (er), aber nicht: Messen des Winkelabstands mit Fotos

• Wenn das hilft, kann ich Pixelpositionen mit mehr Sternen angeben. Ich glaube, 3 sollte ausreichen, aber ich könnte mich irren.

• Ich habe versucht, 3 Sterne auszuwählen, die über das Bild "verteilt" waren (weil ich denke, dass dies den Fehler reduziert, nicht sicher), bin mir aber nicht sicher, ob es mir gelungen ist.

• Ich mache das für mehrere Bilder und möchte eine allgemeine Methode.

• Auf diese Weise kann ich schwächere Sterne / Messier-Objekte / usw. auf dem Bild identifizieren.

• Ich bin sicher, dass viele Astrofotografen dies tun möchten, aber keine vorhandene Software gefunden haben, die dies tut.

EDIT: Danke, whuber! Die gnomonische Projektion hat mir gefehlt. Ich hatte dies bereits unter der Annahme einer linearen Transformation getan:

(* convert RA/DEC to xyz coords on celestial psuedo-sphere of radius 1 *) 
radecxyz[ra_,dec_] = 
{Cos[ra/12*Pi]*Cos[dec/180*Pi],Sin[ra/12*Pi]*Cos[dec/180*Pi],Sin[dec/180*Pi]}; 

(* I no longer have any idea how this works *) 
astrosolve[x_,y_,z_,xwid_,ywid_] := Module[{a,m,ans,nullans}, 
m=Array[a,{2,3}]; 
temp=Solve[{ 
m.radecxyz[x[[1]],x[[2]]]=={x[[3]]-xwid/2,x[[4]]-ywid/2}, 
m.radecxyz[y[[1]],y[[2]]]=={y[[3]]-xwid/2,y[[4]]-ywid/2}, 
m.radecxyz[z[[1]],z[[2]]]=={z[[3]]-xwid/2,z[[4]]-ywid/2} 
}]; 
ans = m /. Flatten[temp]; 
nullans=Flatten[NullSpace[ans]]; 
If[nullans.radecxyz[x[[1]],x[[2]]]<0,nullans=-nullans]; 
Return[{ans,nullans}]; 
]; 

Dabei waren x, y und z jeweils 4-Element-Listen, die aus einem Stern RA, einer Deklination, einer x-Koordinate auf dem Bild und einer y-Koordinate auf dem Bild bestehen. xwid und ywid sind die Breite und Höhe des Bildes. In diesem Fall:

astrosolve[ 
 {19.8463886110, 8.8683219443, 782, 532}, 
 {20.6606352777, 15.9120805555, 311, 146}, 
 {19.4249723610, 3.1147752777, 1023, 815}, 
 1443, 998] 

{ 
 {{-2250.51, -1182.52, 385.689},  {-166.12, -543.746, -2376.73}},  
 {0.480698, -0.861509, 0.163497} 
} 

Unter Bezugnahme auf "{-2250.51, -1182.52, 385.689}" als $ frow, "{-166.12, -543.746, -2376.73}" als $ srow und "{0.480698, -0.861509, 0.163497}" als $ null, Diese PHP-Subroutine übersetzt RA / DEC in xy-Koordinaten:

# radecxy(ra,dec): converts ra/dec to x,y using a quasi-linear transformation 

function radecxy($ra,$dec) { 
    global $null,$frow,$srow,$xwid,$ywid; 
    list($x,$y,$z)=array(cos($dec)*cos($ra),cos($dec)*sin($ra),sin($dec)); 

    $dotprod=$null[0]*$x+$null[1]*$y+$null[2]*$z; 
    if ($dotprod<0) {return(array(-1,-1));}

 list($fx,$fy)  = array($frow[0]*$x+$frow[1]*$y+$frow[2]*$z,$srow[0]*$x+$srow[1]*$y+$srow[2]*$z); 
    $fx+=$xwid/2; 
    $fy+=$ywid/2; 
    if ($fx<0 || $fy<0 || $fx>$xwid || $fy>$ywid) { 
        return(array(-1,-1)); 
    } else { 
        return(array($fx,$fy)); 
    } 
} 

Leider habe ich keine Ahnung mehr, warum dies funktioniert, aber wenn Sie es verwenden + bekannte Sternpositionen hinzufügen, erhalten Sie tolerierbare Ergebnisse (verwenden Sie "Bild anzeigen", um es in voller Größe zu sehen):

Wie Sie jedoch sehen können, sind die Ergebnisse nicht perfekt, was mich davon überzeugt, dass eine lineare Transformation nicht die richtige Antwort war. Ich denke, Gnomonisch könnte der Gral sein, den ich suchte.

Ich werde einen strengen Ansatz skizzieren und angeben, welche Software dabei helfen kann. Das meiste davon wird tangential zu den Interessen der Fotostelle sein, aber da es einige nützliche Erkenntnisse gibt, die für jeden Umstand gelten, unter dem Orte anhand von Messungen an einem Bild geschätzt werden, scheint diese Seite ein vernünftiger Ort für eine solche Analyse zu sein.

Bei der Aufnahme eines Bildes (mit einer Linse, deren Verzerrung korrigiert wurde) wird die Himmelskugel durch den Brennpunkt der Linse auf die Ebene des Sensors projiziert. Dies ist ein schiefer Aspekt einer gnomonischen Projektion .

Mathematisch erfolgt die Konvertierung von (RA, DEC) in mehreren Schritten:

1. Konvertieren Sie (RA, DEC) in sphärische Koordinaten. RA muss von Stunden-Minuten-Sekunden in Grad (oder Bogenmaß) und DEC von Grad-Minuten-Sekunden in Grad (oder Bogenmaß) umgerechnet werden, wobei zu berücksichtigen ist, dass es sich um eine Höhe über dem Flugzeug handelt und nicht um einen Winkel vom Nordpol (Dies ist die übliche sphärische Koordinatenkonvention). Beide Konvertierungen sind einfache Arithmetik.

2. Berechnen Sie (x, y, z) -Koordinaten für die Kugelkoordinaten der Sterne. Dies ist ein StandardStandardkoordinatenumwandlung (einschließlich einfacher Trigonometrie).

3. Drehen Sie die Himmelskugel, um ihre Pole an der Linsenachse auszurichten. Dies ist eine lineare Transformation.

4. Drehen Sie die Himmelskugel um ihre Pole, um sie an die Ausrichtung der Kamera anzupassen (eine weitere lineare Transformation).

5. Platzieren Sie die Bildebene in einer konstanten Höhe z über dem Brennpunkt und ziehen Sie Lichtstrahlen von den Sternen bei (x, y, z) durch den Brennpunkt, bis sie die Ebene abfangen. (Dies ist die gnomonische Projektion und von Natur aus projektiv und nicht linear.)

[In der Figur, die ein planarer Querschnitt durch die Achse der Linse sein soll,

• A ist der Schwerpunkt.
• Halbkreis BCD ist der sichtbare Teil der Himmelskugel.
• Wechselstrompunkte entlang der Linsenachse.
• E, F und G sind Sternpositionen.
• EE, FF und GG sind ihre entsprechenden Positionen auf der (unsichtbaren) Himmelssphäre.
• E ', F' und G 'sind ihre Bilder auf dem Sensor KL (so dass EE', FF 'und GG' Lichtstrahlengänge von den Sternen zum Sensor sind).
• AD ist der Horizont, ab dem die Deklination gemessen wird.
• Alpha ist die Deklination des Sterns E (oder äquivalent eine Winkelkoordinate von EE). Die Sterne F und G haben ähnliche Deklinationen (nicht gezeigt).

Unsere Aufgabe ist es, die mathematische Beziehung zwischen den Winkelkoordinaten für E, F und G - von denen angenommen wird, dass sie mit hoher Genauigkeit bekannt sind - wie Alpha und den gemessenen Koordinaten ihrer Bilder E ', F' und G ' in Pixel entlang des Sensors. Einmal gefunden, kann diese Beziehung wie unten beschrieben invertiert werden, um Winkelkoordinaten von Himmelsobjekten aus den Positionen ihrer Bilder auf dem Sensor abzuschätzen. Der Einfachheit halber ist die Vergrößerung der Linse nicht gezeigt. Bei einer verzerrungsfreien Linse werden die Koordinaten von E ', F' und G 'relativ zur Mitte des Sensors gleichmäßig neu skaliert.]

Dieses Verfahren beschreibt, wie das Licht von einem Stern auf den Sensor gelangt, um eine perfekte einfache Linse zu erhalten. Es handelt sich um diese (unbekannten) Parameter, die bestimmt werden müssen:

• Drei Winkel in (3) und (4) beschreiben die Ausrichtung des Objektivs und der Kamera.

• Ein Skalierungsfaktor in (5), der die kombinierten Effekte der Sensorgröße, der Entfernung vom Brennpunkt und der Vergrößerung der Linse beschreibt.

Aufgrund der Projektion (5) ist dies im Allgemeinen eine komplexe, nichtlineare Transformation, die jedoch eine eindeutige mathematische Beschreibung hat. Wenn wir x = (RA, DEC) die Position eines Sterns bezeichnen lassen, sei Theta die vier Parameter für den Bildgebungsprozess darstellen und y = (Spalte, Zeile) die Pixelkoordinaten darstellen lassen, können wir abstrakt, aber einfacher schreiben

y = f (x, Theta).

Als nächstes - und das ist sehr wichtig - müssen wir Fehler berücksichtigen. Die abgebildeten Sterne befinden sich nicht an genauen Orten. Daher müssen wir einen Fehlerterm in unsere Formel aufnehmen, und es ist üblich (seit etwa 1800), diesen Fehler probabilistisch zu modellieren. Die neue Formel lautet

y = f (x, Theta) + e

Wenn die Linse verzerrungsfrei ist, ist der erwartete Wert von e 0 und seine Standardabweichung ( Sigma ) misst die typische Fehlergröße. Es ist vernünftig anzunehmen, dass die e ungefähr normalverteilt sind, mit ungefähr gleichen Standardabweichungen (was nicht wahr ist, aber für eine erste Analyse ist es eine vernünftige Annahme) und wir können hoffen, dass diese Fehler statistisch unabhängig voneinander sind (was wiederum nicht wahr ist aber es ist eine gute Ausgangsannahme). Dies rechtfertigt eine Lösung der kleinsten Quadrate mit maximaler Wahrscheinlichkeit . Bis zu einer universellen Konstante, deren Wert wir nicht kennen müssen, ist die logarithmische Wahrscheinlichkeit einer bestimmten Beobachtung (x, y) gleich

- | f (x, Theta) - y | ^ 2 / (2 Sigma ^ 2) - 2 log (Sigma).

(Die Absolutwertbalken bezeichnen den euklidischen Abstand in der Abbildungsebene, der wie üblich mit dem Satz von Pythagoras berechnet wird.)

Aufgrund der angenommenen Unabhängigkeit von Fehlern ist die logarithmische Wahrscheinlichkeit des Datensatzes für ein Bild die Summe dieser Werte. Dies ist die "Protokollwahrscheinlichkeit". Die Maximum Likelihood (ML) -Schätzungen der Parameter Theta und Sigma (insgesamt fünf Zahlen) sind diejenigen Werte, die die Log-Wahrscheinlichkeit maximieren.

Wir können und sollten noch weiter gehen. Die Theorie der ML zeigt auch, wie Konfidenzintervalle für die Schätzungen erhalten werden. Intuitiv erzeugen die Fehler in unseren Beobachtungen eine kleine Unsicherheit in den Gelenkwerten der Winkel, des Skalierungsfaktors und der Standardabweichung. Wir benötigen diese Werte, um RA und DEC für alle Pixel in unserem Bild zu schätzen. Durch die Verwendung unsicherer Werte, die unvermeidlich sind, erhalten wir unsichere Ergebnisse. Wenn wir ein Pixel in unserem Bild identifizieren, indem wir einen diffusen Lichtfleck betrachten (der über ungefähr pi * sigma ^ 2 Pixel insgesamt gestreut ist), gibt es außerdem zusätzliche Unsicherheit in den Pixelkoordinaten. Zusammengenommen verbinden sich diese beiden Formen der Unsicherheit. Dies impliziertDie Nettounsicherheit bei der Schätzung von RA und DEC eines Lichtflecks auf dem Bild ist größer als Sie vermuten.

Wenn Sie schließlich eine Messung aus dem Bild vornehmen und damit die wahren Koordinaten eines Sterns oder Himmelsobjekts schätzen, führen Sie eine inverse Regression durch , eine Form der Instrumentenkalibrierung. Inverse Regression ist ein Verfahren, um die soeben beschriebenen Unsicherheiten zu berücksichtigen. Die Ausgabe enthält natürlich die geschätzten Sternkoordinaten für jeden Pixelklecks auf dem Bild. Es enthält auch einen Koordinatenring um diese Schätzung, der auch mit der Position dieses Blobs übereinstimmt. (Dies ist ein gemeinsames "inverses Vorhersageintervall" oder eine Reihe von "Bezugsgrenzen" für RA und DEC des Blobs.) Wenn Sie in der Praxis einen Katalog von Himmelsobjekten konsultieren, können Sie diesen Ring verwenden, um nach allen bekannten Objekten zu suchen, die stimmen mit den Informationen in Ihrem Bild überein. Dies kann eindeutig wertvoller sein als ein vereinfachtes Verfahren, bei dem - manchmal falsch - nur ein einziger Satz von Koordinaten geschätzt wird.

Zusammenfassend wird hier Software benötigt

• Führen Sie die für ML erforderliche nichtlineare Optimierung durch.

• Schätzen Sie Standardfehler in den Schätzungen.

• Führen Sie eine inverse Regression durch.

Kenntnisse mit geeigneter Software, wie dem ML-Befehl von Stata oder Mathematica , sind unerlässlich, wenn Sie dies selbst codieren.

Unabhängig von Ihrem Fachwissen können Sie hier einige Schlussfolgerungen für Ihre Bildgebungsstrategien ziehen:

• Die Genauigkeit des Bildes zum Erhalten einer Fixierung auf einem Objekt kann niemals größer sein als die inhärente Ungenauigkeit in der Bildgebung (gemessen durch Sigma , die typische Größe eines Lichtpunkts auf dem Bild).

• Sie können sich dieser Genauigkeit annähern, indem Sie viele bekannte Sterne identifizieren , nicht nur drei. Dies reduziert die Unsicherheit bei der Umwandlung von Himmel in Bild fast auf Null, wenn sich genügend bekannte Sterne im Bild befinden.

• Es ist richtig, dass die Referenzsterne über das Bild verteilt werden sollen. Es ist auch wichtig, dass sie nicht in einer Reihe stehen (was leider bei den drei in der Frage angegebenen Orten der Fall ist). Wenn Sie es sich leisten können, nur drei Sterne zu finden, bringen Sie sie in ein schönes Dreieck. Wenn sich die Sterne ausrichten, zeigt die statistische Analyse, dass eine große Unsicherheit über Orte in Richtungen senkrecht zur Linie besteht. In diesem speziellen Beispiel ist der geschätzte Fehler ( Sigma ) Hunderte von Pixeln breit. Wenn Sie einen Stern mehr haben, um ein gutes Dreieck zu bilden, sollte dieser Fehler auf ein oder zwei Pixel reduziert werden.

Einige Abschiedsgedanken:

• Durch eine umfassendere statistische Analyse können Linsenaberrationen erkannt und sogar korrigiert werden. Die Idee ist, Abweichungen zwischen erwarteten und tatsächlichen Positionen der Sterne auf dem Bild aufzuzeichnen. Dies ist vergleichbar mit "Warping" - oder " Georeferenzierung " -Kartendaten. Als schnelle und schmutzige Lösung können Sie GIS oder Bildverarbeitungssoftware (wie ENVI ) in Betrieb nehmen, um jedes Bild zu georeferenzieren (oder astroreferenzieren). Eine solche Software führt normalerweise keine ML-Schätzungen projektiver Transformationen durch, kann jedoch Polynomnäherungen höherer Ordnung durchführen. Abhängig von Ihrer Anwendung kann eine Polynomtransformation der Ordnung 2 oder 3 ausreichend funktionieren.

• Es ist möglich, die Genauigkeit zu verbessern, indem mehrere Bilder derselben Objekte kombiniert werden.

Dies mit der gleichen Präzision zu tun, die professionelle Astronomen tun, wäre in der Tat schwierig. Es würde eine äußerst genaue Charakterisierung der von Ihrem Objektiv verursachten Verzerrungen und der Unvollkommenheiten im Sensor Ihrer Kamera erfordern. Allerdings benötigen Sie diesen Genauigkeitsgrad wahrscheinlich nicht. Es sollte für Sie ausreichend sein anzunehmen, dass Ihr Objektiv keine großen Verzerrungen hervorruft (was eine gute Annahme für ein Qualitätsobjektiv ist) und dass Ihr Kamerasensor ziemlich nahe an einem vollkommen regelmäßigen Raster liegt (was eine sehr gute Annahme für ist sogar eine billige Kamera).

Alles, was bleibt, ist die Koordinatentransformation zu erarbeiten, die die Ausrichtung der Kamera beschreibt, dh die Richtung, in die sie gerichtet war, und den Grad, in dem sie gedreht wurde.

Was Sie dann suchen, wird als affine Transformation oder affine Karte bezeichnet. Dies ist nur ein ausgefallener Name für eine Matrix, mit der Sie Ihre Pixelkoordinaten multiplizieren würden, um Ihre astronomischen Koordinaten zu erhalten. Im Fall einer affinen Karte kann diese Transformation einen beliebigen Grad an Rotation, Skalierung, Scherung und Translation umfassen.

Die Bedeutung der Rotationskomponente ist ziemlich offensichtlich. Der Skalierungsfaktor beschreibt einfach, wie viel des Himmels von jedem Pixel in Bezug auf RA / Dez bedeckt ist. Scherung ist eine Transformation, die das Bild eines Rechtecks ​​zu einem Parallelogramm machen würde, aber es sollte keinen dieser Effekte in einem Bild von Objekten im Unendlichen (wie Sternen) geben. Schließlich fügt die Übersetzungskomponente einfach einen Versatz hinzu, um die Tatsache zu berücksichtigen, dass das Pixel (x = 0, y = 0) in Ihrem Bild wahrscheinlich nicht (RA = 0, Dec = 0) entspricht.

Da Ihr Bild 3 Referenzsterne enthält, verfügen Sie über genügend Informationen, um die Beziehung zwischen Ihren Pixelkoordinaten und dem gesuchten RA / Dec zu berechnen. Dies würde durch lineare Anpassung der kleinsten Quadrate (nicht wie oben erwähnt nichtlineare kleinste Quadrate) erfolgen, um die Werte der Matrixkomponenten zu bestimmen, die Ihren Pixelkoordinaten am besten mit der bekannten RA / Dez der Referenzsterne übereinstimmen. Sobald die Matrix erstellt ist, können Sie sie auf die Pixelkoordinaten anderer Sterne anwenden, um deren RA / Dez zu erhalten.

Obwohl ich dies relativ einfach tun könnte, bin ich mir leider nicht sicher, wie ich Ihnen dabei helfen kann. Es würde einige mathematische Fähigkeiten beinhalten, die etwas über den Rahmen von photo.SE hinausgehen. Ich bin Optiker, aber kein großer Fotograf. Die Software, die ich dafür verwenden würde, wurde für Ingenieure entwickelt, um numerische Hochleistungsberechnungen durchzuführen, und ist überhaupt kein fotografisches Werkzeug. Es gibt möglicherweise Möglichkeiten, dies mit Softwarepaketen zu tun, die sich an Fotografen richten, aber ich weiß nichts darüber.