//-----------------------------------------------------------------------------
// NastGrid2d.cpp
//-----------------------------------------------------------------------------

#include "NastGrid2d.h"
#include "NastDebug.h"
#include "NastBox2d.h" 
#include "string.h"

#ifdef _MSC_VER
#   include "strstrea.h"
#else
#   include "strstream.h"
#endif

//-----------------------------------------------------------------------------
//                    Konstruktor + Destruktor
//-----------------------------------------------------------------------------

CNastGrid2d::CNastGrid2d( const CNastBox2d &box,
        		  int nSizeI, int nSizeJ )
    :CNastArray2d<double>( nSizeI, nSizeJ ),
     m_box( box )
{
    NAST_ASSERT( nSizeI >= 2 );
    NAST_ASSERT( nSizeJ >= 2 );

    m_dx = m_box.sizeX() / (sizeI()-1);
    m_dy = m_box.sizeY() / (sizeJ()-1);

    m_rdx = (sizeI()-1) / m_box.sizeX();
    m_rdy = (sizeJ()-1) / m_box.sizeY();

    CNastArray2d<double>::setAll(0.0);
}

CNastGrid2d::CNastGrid2d( const CNastPoint2d &pntMin,
        		  const CNastPoint2d &pntMax,
        		  int nSizeI, int nSizeJ )
    :CNastArray2d<double>( nSizeI, nSizeJ ),
     m_box( pntMin, pntMax )
{
    NAST_ASSERT( nSizeI >= 2 );
    NAST_ASSERT( nSizeJ >= 2 );

    m_dx = m_box.sizeX() / (sizeI()-1);
    m_dy = m_box.sizeY() / (sizeJ()-1);

    m_rdx = (sizeI()-1) / m_box.sizeX();
    m_rdy = (sizeJ()-1) / m_box.sizeY();

    CNastArray2d<double>::setAll(0.0);
}

//  Destruktor
CNastGrid2d::~CNastGrid2d()
{
    NAST_ASSERT_VALID( this );
}

//  Copyconstr.
CNastGrid2d::CNastGrid2d( const CNastGrid2d &other )
    :CNastArray2d<double>( other ),
     m_dx( other.m_dx ),
     m_dy( other.m_dy ),
     m_box( other.m_box )
{}

//  Zuweisungsoperator
const CNastGrid2d&
CNastGrid2d::operator=( const CNastGrid2d &other )
{
    if( this != &other )		//  keine Zuweisung auf sich selbst ?
    {
        m_box     = other.m_box;
        m_dx      = other.m_dx;
        m_dy      = other.m_dy;
    }
    
    return *this;
}

//  absMax()    maximaler Absolutwert
double CNastGrid2d::absMax() const
{
    double max = -1;

    NAST_ASSERT_VALID(this);
    for( int i = 0; i < sizeI(); i++)
        for( int j = 0; j < sizeJ(); j++)
            if( fabs( (*this)(i,j)) > max)
        	max = (*this)(i,j);

    NAST_ASSERT( max >= 0 );	// keine Elemente vorhanden ?

    return max;
}


//  bilineare Interpolation
double CNastGrid2d::interpolate( const CNastPoint2d &pnt ) const
{
    NAST_ASSERT_VALID( this );
    NAST_ASSERT( m_box.isInside( pnt )); //  innerhalb des Gitters 
    
    // lokale Koordinaten des Punktes
    double x = pnt.x() - m_box.pntMin().x();
    double y = pnt.y() - m_box.pntMin().y();

    // Zelle bestimmen
    int i = (int)(floor(x * m_rdx));	 // x / m_dx 
    int j = (int)(floor(y * m_rdy));	 // y / m_dy

    // Koordinaten der Werte
    const double x0 =  i      * m_dx;
    const double x1 = (i + 1) * m_dx;

    const double y0 =  j      * m_dy;
    const double y1 = (j + 1) * m_dy;

    //  Wenn der Punkt auf einer Kante zwischen zwei Werten liegt,
    //  dann gibt es eine Sonderbehandlung. Damit laesst sich auch
    //  gleich der Fall abfangen, dass die der Punkt auf den auesseren
    //  Kanten liegt und somit bei der Interpolation eine Indexueber-
    //  schreitung stattfinden wuerde.

    // Um den Code stabiler zu machen werden kleine Ungenauigkeit 
    // zugelasse
    if( fabs(x0 - x) < NAST_GEOMETRIE_EPS &&    // liegt der Punkt auf einem
        fabs(y0 - y) < NAST_GEOMETRIE_EPS )     // Gitterpunkt ?
    {
        return (*this)(i, j);	// dann nehmen wir den doch als Ergebnis
    }
    else if( fabs(x0 - x) < NAST_GEOMETRIE_EPS)	 // Punkt auf einer vertikalen Kante ?
    {
        return
            ((y1 - y ) * (*this)( i, j   ) +
             (y  - y0) * (*this)( i, j+1 )) * m_rdy;
    }
    else if( fabs(y0 - y) < NAST_GEOMETRIE_EPS ) // Punkt auf einer horizontalen Kante ?
    {
        return 
            ((x1 - x ) * (*this)( i  , j ) +
             (x  - x0) * (*this)( i+1, j )) * m_rdx;
    }

    // ansonsten alle vier Gitterwerte verwenden
    return 
        ((x1 - x ) * (y1 - y ) * (*this)( i  , j   ) +
         (x  - x0) * (y1 - y ) * (*this)( i+1, j   ) +
         (x1 - x ) * (y  - y0) * (*this)( i  , j+1 ) +
         (x  - x0) * (y  - y0) * (*this)( i+1, j+1 )) * m_rdx * m_rdy;
}

// ----------------------------------------------------------------------------
//                                    Debug
// ----------------------------------------------------------------------------


// liefert einen string mit n Leerzeichen
static const char* spaces( unsigned int n )
{
    static const char szSpaces[] = "                                           ";
    NAST_ASSERT( n < strlen( szSpaces ));
    return szSpaces + strlen(szSpaces) - n;
}


// liefert einen string mit n Leerzeichen
static const char* seperate( unsigned int n )
{
    static const char szSeperate[] = "-------------------------------------------";
    NAST_ASSERT( n < strlen( szSeperate ));
    return szSeperate + strlen(szSeperate) - n;
}


// liefert die Zahl n in einem Feld der Groesse nWidth rechtszentriert
// Hmm, ist ein ziemlicher Hack, sollte aber fuer den Zweck reichen
static const char* field( int n , int nWidth )
{
    static char szOut[101];
    ostrstream stream( szOut, 100 );
    stream.width( nWidth);
    stream << n << '\0';
    return stream.str();
 }

//  und noch so eine unschoene Sache
static const char* field( double d, int nWidth )
{
    static char szOut[101];
    szOut[0] = '\0';

    ostrstream stream( szOut, 100 );

    int n = nWidth - 3  - 4; //  -1.234e+01
                             //     ^^^ precision
    NAST_ASSERT( nWidth > 0 );

    if( d == 0)
    {
        stream << spaces(nWidth);
    }
    else if( d >= 0)
    {
        stream.setf( ios::scientific);
        stream.precision(n);
        stream << " " << d << '\0';
    }
    else
    {
        stream.setf( ios::scientific);
        stream.precision(n);
        stream << d << '\0';
    }

    return stream.str();
}

void 
CNastGrid2d::debugDump( CNastDumpContext &dumpContext ) const
{
    CNastObject::debugDump( dumpContext );
    dumpContext << "\tCNastGrid2d" << "\n";
    

    dumpContext << "\t\tUrsprung  x: " << box().minX() << "\n";
    dumpContext << "\t\t          y: " << box().minY() << "\n";

    dumpContext << "\t\tAnzahl    i: " << sizeI() << "\n";
    dumpContext << "\t\t          j: " << sizeJ() << "\n";

    dumpContext << "\t\tGroesse   x: " << m_box.sizeX() << "\n"; 
    dumpContext << "\t\t          y: " << m_box.sizeY() << "\n"; 

    dumpContext << "\t\tAbstand  dx: " << dx() << "\n";
    dumpContext << "\t\t         dy: " << dy() << "\n";

    // Falls das Gitter nicht zu gross ist
    if( sizeI() < 20 && sizeJ() < 20 )
    {
        int i,j;
        const int nSizeBorder  =  4; // Anz. Zeichen in der linken Spalte
        const int nSizeElement =  10; // Anz. Zeichen die ein Eintrag haben soll
        const int nSizeSpace   =  1; // Anz. Zeichen die fuer Trennzeichen verwendet werden
        
        // Zeile mit den Spaltennummern
        dumpContext << spaces( nSizeBorder ) << "|";
        for( i = 0; i < sizeI(); i++)
        {
            dumpContext << field( i, nSizeElement / 2);
            dumpContext << spaces( nSizeElement / 2 + nSizeSpace );
        }
        dumpContext << "\n";
        
        // Zeile mit "----"
        dumpContext << seperate( nSizeBorder );
        for( i = 0; i < sizeI(); i++)
            dumpContext << seperate( nSizeElement + nSizeSpace );
        dumpContext << "\n";
        
        // die Werte
        for( j = sizeJ()-1; j >= 0; j--)
        {
            dumpContext << field( j, nSizeBorder-1) << " |";
            for( i = 0; i < sizeI(); i++)
            {
//          dumpContext << field( operator()(i,j), nSizeElement - nSizeSpace );
        	dumpContext << field( operator()(i,j), nSizeElement );
        	dumpContext << spaces( nSizeSpace );
            }
            dumpContext << "\n";
        }
    }
}

void 
CNastGrid2d::assertValid() const
{
    // zuerst einmal AssertValid der Basisklasse aufrufen
    CNastObject::assertValid(); 
}