=== Arbeitsgruppe 1: Freibord ===

Die Aufgabe der Arbeitsgruppe bestand darin, die I
 
'''Daten''' 

Der Gruppe stand ein ASAR-Satellitenbild zur Verfügung, das einen Auschnitt des Weddellmeeres zeigt und ein ICESat-Datensatz, der die geographischen Positionen (lon, lat) und die an diesen Punkten gemessenen Freibordhöhen (cm) für einen Überflug quer durch den ASAR-Ausschnitt beinhaltet.

{{attachment:schema1.jpg}}

 
'''Methodik''' 
 
Theorie zur Koordinatentransformation:

Eine affine Abbildung ist eine lineare Koordinatentransformation, die die elementaren Transformationen Translation, Rotation, Dilatation, Stauchung und Scherung umfasst. Sie kann durch Vektoraddition und Matrixmultiplikation ausgedrueckt werden:
{{{#!latex
\[\left(\begin{array}{c} x'  \\ y'  \end{array}\right) =  
\left(\begin{array}[c]{ccc} a_{\rm{11}} & a_{\rm{12}}\\a_{\rm{21}} & a_{\rm{22}}\end{array}\right)
\left(\begin{array}{c} x  \\ y  \end{array}\right)
+
\left(\begin{array}{c} t_x  \\ t_y  \end{array}\right)\]
}}}
Homogene Koordinaten:
{{{#!latex
\[\left(\begin{array}{c} x'  \\ y' \\ 1 \end{array}\right) =  
\left(\begin{array}[c]{ccc} a_{\rm{11}} & a_{\rm{12}} & t_x \\ a_{\rm{21}} & a_{\rm{22}} & t_y \\ 0 & 0 & 1 \end{array}\right)
\left(\begin{array}{c} x  \\ y \\ 1 \end{array}\right)\]
}}}
Für drei nichtkollineare Punkte ergibt sich damit folgendes Gleichungssystem:
{{{#!latex
\[\left(\begin{array}[c]{ccc} x_1' & x_2' & x_3' \\ y_1' & y_2' & y_3' \\ 1 & 1 & 1 \end{array}\right)=
\left(\begin{array}[c]{ccc} a_{\rm{11}} & a_{\rm{12}} & t_x \\ a_{\rm{21}} & a_{\rm{22}} & t_y \\ 0 & 0 & 1 \end{array}\right)
\left(\begin{array}[c]{ccc} x_1 & x_2 & x_3 \\ y_1 & y_2 & y_3 \\ 1 & 1 & 1 \end{array}\right)\]
}}}
oder
{{{#!latex 
\[P'=AP\]    
}}}
Die Transformationsmatrix A für drei nichtkollineare Punkte lässt sich dann einfach aus 
{{{#!latex 
\[A=P'P^{-1}\]
}}}
bestimmen.

Transformationskoeffizienten für mehr als drei nichtkollineare Punkte erhält man mit der Methode der kleinsten Quadrate aus:
{{{#!latex 
\[A=P'P^T(PP^T)^{-1}\]
}}}
 
(für mehr Informationen siehe B. Jähne, Digitale Bildverarbeitung, Kapitel 10.4)


Arbeitsschritte:

Die Gruppe hat Programme bzw. Funktionen erarbeitet, die folgendes tun:
 * Einlesen der ICESat-Datei und umrechnen der geographischen Koordinaten des Ueberfluges in polarstereographische Koordinaten.
 * Die Eckpunkte des ASAR-Bildes, die von der Arbeitsgruppe 0 in geographischen Koordinaten übergeben wurden, werden ebenfalls in polarstereographische Koordinaten umgerechnet 
 * Die Koordinatentransformation wird so durchgeführt, dass man als Ergebnis den ICESat-Datensatz als normierte Bildkoordinaten erhält. Dazu werden zunächst die vier ASAR-Eckpunkte in normierte Bildkoordinaten gebracht und anschließend die ICESat-Daten auf dasselbe Koordiantensystem transformiert.

{{attachment:schema.jpg}}
 * Als Endergebnis wird eine Matrix erzeugt, die die Messpositionen des ICESat-Ueberfluges im ASAR-Auschnitt in Bildkoordinaten und die zugeoerigen Freibordoehen enthält.
  
'''Ergebnisse''' 

Der Output besteht dann wie oben gesagt aus einer Matrix mit den Messpositionen des ICESat-Ueberfluges und der zugehoerigen Freibordhoehen.
 
(Output, Statistik)
  
'''Diskussion''' 


=== Programme ===

'''fbh_bildkoordinaten.py'''

{{{#!python
from polar_projection import *
from read_asar import *
from read_icesat import *
from coord_transform import *
from scipy import *

def fit_freeboard_ASAR(filename1,filename2):
    """filename1: ASAR data file, filename2: freeboard data file
       creates new coordinate system defined by corners of ASAR image and selects freeboard values within ASAR image box
       returns an array containing normalized image coordinates and corresponding freeboard values:
       [x_coordinate, y_coordinate, freeboardheight(cm)]"""
    
    sgn=-1  #Antarctica
    ASAR=array(read_asar_corners(filename1))
    ASAR_p=zeros(8)
    for k in arange(0,7,2):       #computing polarstereographic coordinates
        ASAR_p[k:k+2]=mapll(ASAR[k],ASAR[k+1],sgn)

    A=coord_transformation(ASAR_p) #computing transformation matrix A for coordinate
                                   #transformation into image coordinates
 
    # reading freeboard data and computing geographic into polarstereographic coordinates
    ICESAT_p,fbh=read_icesat(filename2,sgn)    #fbh are measured freeboard heights in cm

    # calculating new coordinates for freeboard data
    x_neu=[]
    y_neu=[]
    for x,y in zip(ICESAT_p[0],ICESAT_p[1]):
        x_neu.append(dot(array([A[0,0],A[0,1]]),array([x,y]))+A[0,2])
        y_neu.append(dot(array([A[1,0],A[1,1]]),array([x,y]))+A[1,2])
    x_n=array(x_neu)
    y_n=array(y_neu)
    fbh_n=array(fbh)

    x_n_limited=clip(x_n,0.,1.)
    x_indices=nonzero(x_n==x_n_limited)
    x_xind=x_n[x_indices]
    y_xind=y_n[x_indices]
    fbh_xind=fbh_n[x_indices]
    
    y_n_limited=clip(y_xind,0.,1.)
    y_indices=nonzero(y_xind==y_n_limited)
    x_bild=x_xind[y_indices]
    y_bild=y_xind[y_indices]
    fbh_bild=fbh_xind[y_indices]
    
    x_y_fbh=array([x_bild,y_bild,fbh_bild])

    return x_y_fbh
}}}

'''coord_transform.py'''
{{{#!python
from scipy import linalg as la

def coord_transformation(ASAR_p):
    # polarstereographic coordinate system
    y00,x00,y01,x01,y02,x02,y03,x03=int(ASAR_p[1]),int(ASAR_p[0]),int(ASAR_p[7]),int(ASAR_p[6]),int(ASAR_p[5]),int(ASAR_p[4]),int(ASAR_p[3]),int(ASAR_p[2])
    # new coordinate system with normalized coordinates  
    y10,x10,y11,x11,y12,x12,y13,x13=0,0,1,0,1,1,0,1

    # calculating transformation matrix:
    P0=array([[x00, x01, x02, x03],[y00,y01,y02,y03],[1.0,1.0,1.0,1.0]])
    P1=array([[x10, x11, x12, x13],[y10,y11,y12,y13],[1.0,1.0,1.0,1.0]])

    Faktor1=dot(P1,la.transpose(P0))
    Faktor2=la.inverse(dot(P0,la.transpose(P0)))
    A=dot(Faktor1,Faktor2)  # Transformation matrix
    return A
  }}}

'''read_icesat.py'''

{{{#!python
# reading freeboard data

import string
from geo_polar import *
from scipy import io

def read_icesat(filename,sgn):
    data=io.read_array(filename)
    polar=mapll(data[:,1],data[:,0],sgn)
    fbh=data[:,2]
    return polar,fbh
}}}

Die benötigten Module polar_projection.py und read_asar.py sind auf der Seite der Arbeitsgruppe 0 [[AG0_ASAR_Einlesen]] zu finden. 

'''fbh_bildkoordinaten_test.py'''

Zum Testen hängt man an das obige Programm folgende Zeilen an: 
{{{#!python
filename1='ASA_IMP_1PNDPA20060617_043346_000000162048_00362_22460_2136.N1'
filename2='LonLatFre_1706_6.xyz'
ergebnis=fit_freeboard_ASAR(filename1,filename2)
}}}



{{attachment.schemabild2.jpg}}