Subversion Repositories public iLand

#include "hemigrid.h"

//#include <algorithm>

#include <QtCore>

// needed for drawing: image & color

#include <QImage>

#include <QColor>

//////////////////////////////////////////////////////

// Setup memory for the Grid.

// @param size of a single pixel in degree

//////////////////////////////////////////////////////

void HemiGrid::setup(double cellsize_degree)

{

      // setup grid...

      mMatrixCountAzimuth = int(360 / cellsize_degree);

      mMatrixCountElevation = int(90 / cellsize_degree);

      mCellSizeDegree = cellsize_degree;

      // size occupied by one pixel in rad

      mMatrixCellSize = cellsize_degree * M_PI / 180.;

      if (mMatrix)

        delete[] mMatrix;

      mMatrix = new double[mMatrixCountAzimuth*mMatrixCountElevation];

      clear();

}

//////////////////////////////////////////////////////

// Clear the Grid.

// @param SetWith value used to fill (default 0.)

//////////////////////////////////////////////////////

void HemiGrid::clear(double SetWith)

{

    if (mMatrix) {

      // blank matrix...

      for (int i=0;i<mMatrixCountAzimuth*mMatrixCountElevation;i++)

         mMatrix[i]=SetWith;

    }

}

void HemiGrid::matrixMinMax(double &rMatrixMin, double &rMatrixMax) const

{

    rMatrixMin = 100000000.;

    rMatrixMax = -1000000000;

    if (mMatrix) {

      // blank matrix...

      for (int i=0;i<mMatrixCountAzimuth*mMatrixCountElevation;i++) {

         rMatrixMin = std::min(mMatrix[i], rMatrixMin);

         rMatrixMax = std::max(mMatrix[i], rMatrixMax);

      }

    }

}

double HemiGrid::sum(const double elevation) const

{

    int ie=indexElevation(elevation);

    int emax=matrixCountElevation();

    int amax=matrixCountAzimuth();

    double value = 0;

    for (int e=ie;e<emax;e++)

        for (int a=0;a<amax;a++)

            value+=rGetByIndex(a,e);

    return value;

}

//////////////////////////////////////////////////////

// Dump the grid into a TStringList

//////////////////////////////////////////////////////

//void HemiGrid::DumpGrid(TStrings* List)

//{

//    AnsiString Line;

//    for (int i=0;i<mMatrixCountAzimuth;i++)  {

//       Line=AnsiString(i);

//       for (int j=0;j<mMatrixCountElevation;j++)  {

//           Line+=";"+AnsiString(rGetByIndex(i,j));

//       }

//       List->Add(Line);

//    }

//

//}

//////////////////////////////////////////////////////

// Get Gridvalue by integer internal grid indices

// @param iAzimuth index of azimuth (0..count(pixels)-1

// @param iElevation index of elevatin (0..count-1)

// @return ref. to grid-value

//////////////////////////////////////////////////////

double& HemiGrid::rGetByIndex(const int iAzimuth, const int iElevation) const

{

if (iAzimuth < mMatrixCountAzimuth && iElevation < mMatrixCountElevation

        && iAzimuth>=0 && iElevation>=0)

       return mMatrix[iAzimuth*mMatrixCountElevation + iElevation];

    else

       throw QString("TSolar::Rad::Get() - invalid indices");

}

//////////////////////////////////////////////////////

// Get Gridvalue

// @param Azimuth azimuth angle (-pi .. +pi, 0=south)

// @param Elevation elevation angle (0=horizon, pi/2: zenith)

// @return ref. to grid-value

//////////////////////////////////////////////////////

double& HemiGrid::rGet(const double Azimuth, const double Elevation) const

{

    // Azimuth goes from -pi .. +pi -> move to 0..2pi, scale to 0..1 and convert to integer indices

    int iAzimuth = indexAzimuth(Azimuth);

    // Elevation goes from 0..90° = 0..pi/2

    int iElevation = indexElevation(Elevation);

    return rGetByIndex(iAzimuth, iElevation);

}

void HemiGrid::modify(const HemiGrid &Source, const ModifyMode mode)

{

        for (int i=0; i<mMatrixCountAzimuth*mMatrixCountElevation; i++) {

            switch (mode) {

                case Add: mMatrix[i]+=Source.mMatrix[i]; break;

                case Multiply: mMatrix[i] *= Source.mMatrix[i]; break;

                case SetTo: mMatrix[i] = Source.mMatrix[i]; break;

            }

        }

}

/** retrieve total sum of the hemigrid.

  @param Weighter if available (non null) each cell value of this grid is multiplied with the weighter grid (note: no additional calculations are performed) */

double HemiGrid::getSum(const HemiGrid *Weighter) const

{

    double Sum=0.;

    if (Weighter) {

        if (Weighter->matrixCountAzimuth()!=this->matrixCountAzimuth()

            || Weighter->matrixCountElevation() != this->matrixCountElevation())

                throw QString("HemiGrid::getSum: invalid weighing object!");

        for (int i=0; i<mMatrixCountAzimuth*mMatrixCountElevation; i++) {

            Sum += mMatrix[i]*Weighter->mMatrix[i];

        }

        return Sum;

    }

    for (int i=0; i<mMatrixCountAzimuth*mMatrixCountElevation; i++) {

        Sum += mMatrix[i];

    }

    return Sum;

}

void HemiGrid::projectLine(const double &x, const double &y, const double &deltah, const double &r, double &elevation, double &azimuth1, double &azimuth2)

{

   // transform coordinates....

   // distance to footing point (x/y)

   double distance = sqrt(x*x + y*y);

   elevation = atan2(deltah, distance);

   double azimuth = atan2(x,y);

   // distance to point (x/y/z)

   double dist3 = sqrt(x*x+y*y+deltah*deltah);

   double azimuth_delta = atan2(r, dist3);

   azimuth1 = azimuth - azimuth_delta;

   azimuth2 = azimuth + azimuth_delta;

}

//////////////////////////////////////////////////////

// Modify a rect defined by coordinates, a mode and a value

// @param  elow, ehigh: elevation angles

// @param alow1, alow2: azimuth angles at base, alow1 must be < alow2

// @param ahigh1, ahigh2: azimuth angles at the top, ahigh1 must be < ahigh2

// @param ModifyMode multiply or add to the matrix

// @param Value value used for multiply/add,...

//////////////////////////////////////////////////////

void HemiGrid::modifyAngleRect( const double &elow, const double &alow1, const double &alow2,

                    const double &ehigh, const double &ahigh1, const double &ahigh2,

                    const ModifyMode mode, const double &value)

{

    int i_e_low = indexElevation(elow);

    int i_e_high = indexElevation(ehigh);

    int i_a_low1 = indexAzimuth(alow1);

    int i_a_low2 = indexAzimuth(alow2);

    int i_a_high1 = indexAzimuth(ahigh1);

    int i_a_high2 = indexAzimuth(ahigh2);

    int i_a_min = std::min(i_a_low1, i_a_high1);

    int i_a_max = std::max(i_a_low2, i_a_high2);

    for (int e = i_e_low; e<=i_e_high; e++) {

       for (int a = i_a_min; a<=i_a_max; a++) {

          double &ref = rGetByIndex((a+mMatrixCountAzimuth)%mMatrixCountAzimuth, e);

          switch (mode) {

            case Multiply: ref*=value; break;

            case Add: ref+=value; break;

            case SetTo: ref=value; break;

          }

       }

    }

}

//////////////////////////////////////////////////////

// Project a cylinder onto the grid.

// @param deltax, deltay: relative position

// @param offsetz: relative height of lower edge

// @param height, r: height and radius of cylinder

// @param mode, value: mode and value of change

// @param  elow, ehigh: elevation angles

//////////////////////////////////////////////////////

void HemiGrid::projectCylinder(const double &deltax, const double &deltay,

                     const double &offsetz, const double &height, const double &r,

                     const ModifyMode mode, const double &value)

{

    double elow, ehigh; // elevations...

    double al1, al2, ah1, ah2;

    projectLine(deltax, deltay, std::max(offsetz+height, 0.), r, ehigh, ah1, ah2);

    // ignore if top is t

    if (ehigh < 5.*M_PI/180.)

       return;

    projectLine(deltax, deltay, std::max(offsetz, 0.), r, elow, al1, al2);

    modifyAngleRect(elow, al1, al2, ehigh, ah1, ah2, mode, value);

}

void HemiGrid::paintGrid(QImage &image) const

{

    int dx = image.width();

    int dy = image.height();

    int ix,iy;

    double x,y;

    int maxsize = std::min(dx,dy);

    double mmin=0, mmax=0; // min and max values of matrix

    //if (mmax==0) mmax=1;

    matrixMinMax(mmin, mmax);

    QColor col;

    double phi, r, elevation, value;

    for (ix=0;ix<maxsize;ix++){

        for (iy=0;iy<maxsize;iy++){

            // to -1..1

            x = 2*(ix-maxsize/2.) / maxsize;

            y = 2*(iy-maxsize/2.) / maxsize;

            // get phi,r (polar coords)

            if (x==0 && y==0)

                continue;

            r = sqrt(x*x + y*y);

            if (r>=1)

                continue;

            phi = atan2(y,x);

            elevation = M_PI_2 - (r * M_PI_2);

            value = rGet(phi, elevation);

            value /= mmax;

            col = QColor::fromHsvF((1.-value)*0.666666 ,0.9, 0.9); // hue from 0..240 = red to blue

            image.setPixel(ix,maxsize-iy,col.rgb());

        }

    }

}

QString HemiGrid::dumpGrid() const

{

    QString text;

    text="rows: azimuth steps, cols: elevation. (first value in row: index of azimuth)\n";

    for (int i=0;i<mMatrixCountAzimuth;i++)  {

       text+=QString::number(i);

       for (int j=0;j<mMatrixCountElevation;j++)  {

           text+=";"+QString::number(rGetByIndex(i,j));

       }

       text+="\n";

    }

    return text;

}