rectangle.c

Go to the documentation of this file.

#include "rectangle.h"

#include "random.h"
#include "definitions.h"

#include <string.h>
#include <stdlib.h>
#include <math.h>
#include <assert.h>


void InitRectangle(unsigned int dim,double *l,double *u,Trectangle *b)
{
  assert(dim>0);

  b->n=dim;

  NEWZ(b->l,2*b->n,double);
  if (l!=NULL)
    memcpy(b->l,l,b->n*sizeof(double));

  b->u=&(b->l[b->n]);
  if (u!=NULL)
    memcpy(b->u,u,b->n*sizeof(double));
}

void InitRectangleFromPoint(unsigned int dim,double *p,Trectangle *b)
{
  assert(dim>0);

  b->n=dim;

  NEW(b->l,2*b->n,double);
  b->u=&(b->l[b->n]);
  memcpy(b->l,p,b->n*sizeof(double));
  memcpy(b->u,p,b->n*sizeof(double));
}

inline unsigned int GetRectangleDim(Trectangle *b)
{
  return(b->n);
}

void ExpandRectangle(double *p,Trectangle *b)
{
   unsigned int i; 
   
   for(i=0;i<b->n;i++)
     {
       if (p[i]<b->l[i])
         b->l[i]=p[i];
       else
         {
           if (p[i]>b->u[i])
             b->u[i]=p[i];
         }
     }
}

double EnlargeRectangleWithLimits(double r,Trectangle *limits,
                                  Trectangle *bIn,Trectangle *bOut)
{
  unsigned int i; 
  double s,v;

  assert(bIn->n==limits->n);
  assert(bOut->n==bIn->n);

  v=1.0;
  for(i=0;i<bIn->n;i++)
    {
      s=bIn->l[i]-r;
      bOut->l[i]=(s>limits->l[i]?s:limits->l[i]);

      s=bIn->u[i]+r;
      bOut->u[i]=(s<limits->u[i]?s:limits->u[i]);

      v*=(bOut->u[i]-bOut->l[i]);
    }
  return(v);
}

void CopyRectangle(Trectangle *b_out,Trectangle *b_in)
{
  b_out->n=b_in->n;

  NEW(b_out->l,2*b_out->n,double);
  memcpy(b_out->l,b_in->l,2*b_out->n*sizeof(double));

  b_out->u=&(b_out->l[b_out->n]);
}

inline double GetRectangleLowerLimit(unsigned int i,Trectangle *b)
{
  assert(i<b->n);

  return(b->l[i]);
}

inline double GetRectangleUpperLimit(unsigned int i,Trectangle *b)
{
  assert(i<b->n);

  return(b->u[i]);
}

inline void GetRectangleLimits(unsigned int i,double *l,double *u,Trectangle *b)
{
  assert(i<b->n);

  *l=b->l[i];
  *u=b->u[i];
}

inline void SetRectangleLowerLimit(unsigned int i,double l,Trectangle *b)
{
  assert(i<b->n);

  b->l[i]=l;
}

inline void SetRectangleUpperLimit(unsigned int i,double u,Trectangle *b)
{
  assert(i<b->n);

  b->u[i]=u;
}

inline void SetRectangleLimits(unsigned int i,double l,double u,Trectangle *b)
{
  assert(i<b->n);

  b->l[i]=l;
  b->u[i]=u;
}

void RandomPointInRectangle(double *c,Trectangle *b)
{
  unsigned int i;
  
  for(i=0;i<b->n;i++)
    c[i]=b->l[i]+randomDouble()*(b->u[i]-b->l[i]);
}

double SquaredDistanceToRectangleDimension(unsigned int dim,double p,unsigned int *tp,Trectangle *b)
{
  double l,u,d,c,d1,d2,r;
  
  assert(dim<b->n);

  l=b->l[dim];
  u=b->u[dim];
  r=(u-l);
  d=0.0;
  if ((tp==NULL)||(tp[dim]==TOPOLOGY_R)||(r>M_2PI+ANGLE_ACCURACY))
    {
      /* topology real, or range larger than 2pi */
      if (p<l)
        d=l-p;
      else
        {
          if (p>u)
            d=p-u;
        }
    }
  else
    {
      if (r<M_2PI-ANGLE_ACCURACY)
        {
          /* range < 2pi */
          PI2PI(l);
          PI2PI(u);
        }
      else
        {
          /* range == 2pi */
          l=-M_PI;
          u=+M_PI;
        }

      c=p;
      PI2PI(c);

      if (u<l)
        {
          /* The interval of interest is [-pi,u] U [l,pi] */
          if ((u<c)&&(c<l))
            {
              d1=c-u;
              d2=l-c;
              d=(d1<d2?d1:d2);
            }
        }
      else
        {
          /* The interval of interest is [l,u] */
          if (c<l)
            {
              d1=l-c;
              if (d1>M_PI) d1=M_2PI-d1;
              d2=u-c;
              if (d2>M_PI) d2=M_2PI-d2;
              d=(d1<d2?d1:d2);
            }
          else
            {
              if (c>u)
                {
                  d1=c-l;
                  if (d1>M_PI) d1=M_2PI-d1;
                  d2=c-u;
                  if (d2>M_PI) d2=M_2PI-d2;
                  d=(d1<d2?d1:d2);
                }
            }
        }
    }
  return(d*d);
}


double SquaredDistanceToRectangle(double t2,double *p,unsigned int *tp,Trectangle *b)
{
  unsigned int i;
  double l,u,d,v,c,d1,d2,r;
  
  d=0.0;
  if (tp==NULL)
    {
      /* topology real */
      for(i=0;((i<b->n)&&(d<t2));i++)
        {
          l=b->l[i];
          if (p[i]<l)
            {
              v=l-p[i];
              d+=v*v;
            }
          else
            {
              u=b->u[i];
              if (p[i]>u)
                {
                  v=p[i]-u;
                  d+=v*v;
                }
            }
        }
    }
  else
    {
      for(i=0;((i<b->n)&&(d<t2));i++)
        {
          l=b->l[i];
          u=b->u[i];

          r=(u-l);
          if ((tp[i]==TOPOLOGY_R)||(r>M_2PI+ANGLE_ACCURACY))
            {
              /* topology real, or range larger than 2pi */
              if (p[i]<l)
                {
                  v=l-p[i];
                  d+=v*v;
                }
              else
                {
                  if (p[i]>u)
                    {
                      v=p[i]-u;
                      d+=v*v;
                    }
                }
            }
          else
            {
              if (r<M_2PI-ANGLE_ACCURACY)
                {
                  /* range < 2pi */
                  PI2PI(l);
                  PI2PI(u);
                }
              else
                {
                  /* range == 2pi */
                  l=-M_PI;
                  u=+M_PI;
                }

              c=p[i];
              PI2PI(c);

              if (u<l)
                {
                  /* The interval of interest is [-pi,u] U [l,pi] */
                  if ((u<c)&&(c<l))
                    {
                      d1=c-u;
                      d2=l-c;
                      v=(d1<d2?d1:d2);
                      d+=v*v;
                    }
                }
              else
                {
                  /* The interval of interest is [l,u] */
                  if (c<l)
                    {
                      d1=l-c;
                      if (d1>M_PI) d1=M_2PI-d1;
                      d2=u-c;
                      if (d2>M_PI) d2=M_2PI-d2;
                      v=(d1<d2?d1:d2);
                      d+=v*v;
                    }
                  else
                    {
                      if (c>u)
                        {
                          d1=c-l;
                          if (d1>M_PI) d1=M_2PI-d1;
                          d2=c-u;
                          if (d2>M_PI) d2=M_2PI-d2;
                          v=(d1<d2?d1:d2);
                          d+=v*v;
                        }
                    }
                }
            }
        }
    }

  return(d);
}

unsigned int GetRectangleSplitDim(Trectangle *b)
{
  unsigned int i;
  unsigned int i_max;
  double s;
  double s_max;
      
  i_max=0;
  s_max=-1.0;
  for(i=0;i<b->n;i++)
    {
      s=(b->u[i]-b->l[i]);
      if (s>s_max)
        {
          s_max=s;
          i_max=i;
        }
    }
  return(i_max);
}

void PrintRectangle(FILE *f,Trectangle *b)
{
  unsigned int i;
  
  fprintf(f,"{ ");
  for(i=0;i<b->n;i++)
    fprintf(f,"[%f,%f] ",b->l[i],b->u[i]);
  fprintf(f,"}");
}

inline void DeleteRectangle(Trectangle *b)
{
  free(b->l);
}