/Users/craigcornelius/Projects/SPRING Mac Release 0.2/triangles.h

Go to the documentation of this file.

// $Id: triangles.h,v 1.18 2002/01/02 01:40:46 kevin Exp $
// $Copyright: (c)2001 National Biocomputation Center, Stanford University $

#ifndef _TRIANGLES_H_
#define _TRIANGLES_H_

#include <stdio.h>

#include "point3d.h"
#include "reallocablearray.h"
#include "space.h"
#include "socket.h"

#define COEFF 0.5
#define WIDTH 1.0

class triangle {
public :
    Point3D p1;
    Point3D p2;
    Point3D p3;
    Point3D t1;
    Point3D t2;
    Point3D t3;
    Point3D normal;
    double d;   // The 4th useful parameters for the equation of the plane.
    char rigidity;
    int force;   // This is true if I had a force last time I ask for getForce.
    
    triangle(const Point3D& np1, const Point3D& np2,
        const Point3D& np3, const Point3D& nt1,
        const Point3D& nt2, const Point3D& nt3,
        char nrigidity) : p1(np1), p2(np2), p3(np3), t1(nt1), t2(nt2), t3(nt3)
           {
        force = 0;
        Point3D u = p2 - p1, v = p3 - p1;
        
        rigidity = nrigidity;
        normal = Point3D(u.y*v.z - u.z*v.y, u.z*v.x - u.x*v.z, u.x*v.y - u.y*v.x);
        normal /= sqrt(normal.Squared());
        d = - (normal.x * p1.x + normal.y * p1.y + normal.z * p1.z);
        if(t1.Dot(normal) < 0)
            t1 = -t1;
        if(t2.Dot(normal) < 0)
            t2 = -t2;
        if(t3.Dot(normal) < 0)
            t3 = -t3;
           }
    triangle() : p1(), p2(), p3(), t1(), t2(), t3(), normal()
    {
        force = 0;
        rigidity = 0;
        d = 0;
    }
    ~triangle()
    {
    }
    
    long gethash()
    {
        long res;
        
        res = hashdouble(p1.x); res ^= hashdouble(p1.y); res ^= hashdouble(p1.z);
        res ^= hashdouble(p2.x); res ^= hashdouble(p2.y); res ^= hashdouble(p2.z);
        res ^= hashdouble(p3.x); res ^= hashdouble(p3.y); res ^= hashdouble(p3.z);
        res ^= hashdouble(t1.x); res ^= hashdouble(t1.y); res ^= hashdouble(t1.z);
        res ^= hashdouble(t2.x); res ^= hashdouble(t2.y); res ^= hashdouble(t2.z);
        res ^= hashdouble(t3.x); res ^= hashdouble(t3.y); res ^= hashdouble(t3.z);
        return res;
    }
    
    Point3D *get_force(Point3D *pos, Point3D *res)
    {
        int test[3];
        Point3D proj;
        double d1,d2,d3;
        double l;
        double sign;
        
        if(res == NULL)
            res = new Point3D();
        
        l = ((normal.Dot(*pos) + d) / normal.Squared());
        
        if(l<0)
        {
            sign = -1;
            l = -l;
        }
        else
            sign = 1;
        
        if(l >= WIDTH)
        {
            force = 0;
            *res = Point3D(0,0,0);
            return res;
        }
        else
        {
            //   cerr << fabs(l) << endl;
        }  
        
        proj = *pos - sign * l * normal;
        
        test[0] = ((p2 - p1).Cross(proj - p1).Dot(normal) > 0?1:-1);
        test[1] = ((p3 - p2).Cross(proj - p2).Dot(normal) > 0?1:-1);
        test[2] = ((p1 - p3).Cross(proj - p3).Dot(normal) > 0?1:-1);
        
        if(test[0] == 0)
            *res = Point3D(t1);
        else if(test[1] == 0)
            *res = Point3D(t2);
        else if(test[2] == 0)
            *res = Point3D(t3);
        else if((test[0] == test[1]) && (test[1] == test[2]))
        {
            d1 = sqrt((proj - p1).Squared());
            d2 = sqrt((proj - p2).Squared());
            d3 = sqrt((proj - p3).Squared());
            *res = Point3D(((1.0/d1) * t1) + ((1.0/d2) * t2) + ((1.0/d3) * t3));
            *res *= 1.0 / ((1.0/d1) + (1.0/d2) + (1.0/d3));
        }
        else
        {
            *res = Point3D(0,0,0);
            force = 0;
            return res;
        }
        force = 1;
        *res *= (1.0 - pow((l/WIDTH), 2.0)) * rigidity * COEFF;
        cerr << "There is a force of magnitude: " << res->Squared() << "." << endl;
        return res;
    }
    
    void print_datas()
    {
    fprintf(stderr, "Triangle class.\nRigidity: %d\np1: %f,%f,%f\nt1: %f,%f,%f\n\
                                                                                  p2: %f,%f,%f\nt2: %f,%f,%f\n\
p3: %f,%f,%f\nt3: %f,%f,%f\n",
rigidity,p1.x,p1.y,p1.z,t1.x,t1.y,t1.z,
p2.x,p2.y,p2.z,t2.x,t2.y,t2.z,
p3.x,p3.y,p3.z,t3.x,t3.y,t3.z);
    }
    
    void receive_data(char *buf)
    {
        int j;
        Point3D u,v;
        
        j=0;
        rigidity = buf[j];
        j++;
        receive_float_from_buf(buf + j, &(p1.x));
        j+=4;
        receive_float_from_buf(buf + j, &(p1.y));
        j+=4;
        receive_float_from_buf(buf + j, &(p1.z));
        j+=4;
        receive_float_from_buf(buf + j, &(t1.x));
        j+=4;
        receive_float_from_buf(buf + j, &(t1.y));
        j+=4;
        receive_float_from_buf(buf + j, &(t1.z));
        j+=4;
        receive_float_from_buf(buf + j, &(p2.x));
        j+=4;
        receive_float_from_buf(buf + j, &(p2.y));
        j+=4;
        receive_float_from_buf(buf + j, &(p2.z));
        j+=4;
        receive_float_from_buf(buf + j, &(t2.x));
        j+=4;
        receive_float_from_buf(buf + j, &(t2.y));
        j+=4;
        receive_float_from_buf(buf + j, &(t2.z));
        j+=4;
        receive_float_from_buf(buf + j, &(p3.x));
        j+=4;
        receive_float_from_buf(buf + j, &(p3.y));
        j+=4;
        receive_float_from_buf(buf + j, &(p3.z));
        j+=4;
        receive_float_from_buf(buf + j, &(t3.x));
        j+=4;
        receive_float_from_buf(buf + j, &(t3.y));
        j+=4;
        receive_float_from_buf(buf + j, &(t3.z));
        u = p2 - p1;
        v = p3 - p1;
        normal = Point3D(u.y*v.z - u.z*v.y, u.z*v.x - u.x*v.z, u.x*v.y - u.y*v.x);
        normal /= sqrt(normal.Squared());
        d = - (normal.x * p1.x + normal.y * p1.y + normal.z * p1.z);
        if(t1.Dot(normal) < 0)
            t1 = -t1;
        if(t2.Dot(normal) < 0)
            t2 = -t2;
        if(t3.Dot(normal) < 0)
            t3 = -t3;
    }
    
    void data_to_send(char *buf)
    {
        int j;
        
        j = 0;
        buf[j] = rigidity;
        j++;
        copy_float_to_buf(buf + j, p1.x);
        j+=4;
        copy_float_to_buf(buf + j, p1.y);
        j+=4;
        copy_float_to_buf(buf + j, p1.z);
        j+=4;
        copy_float_to_buf(buf + j, t1.x);
        j+=4;
        copy_float_to_buf(buf + j, t1.y);
        j+=4;
        copy_float_to_buf(buf + j, t1.z);
        j+=4;
        copy_float_to_buf(buf + j, p2.x);
        j+=4;
        copy_float_to_buf(buf + j, p2.y);
        j+=4;
        copy_float_to_buf(buf + j, p2.z);
        j+=4;
        copy_float_to_buf(buf + j, t2.x);
        j+=4;
        copy_float_to_buf(buf + j, t2.y);
        j+=4;
        copy_float_to_buf(buf + j, t2.z);
        j+=4;
        copy_float_to_buf(buf + j, p3.x);
        j+=4;
        copy_float_to_buf(buf + j, p3.y);
        j+=4;
        copy_float_to_buf(buf + j, p3.z);
        j+=4;
        copy_float_to_buf(buf + j, t3.x);
        j+=4;
        copy_float_to_buf(buf + j, t3.y);
        j+=4;
        copy_float_to_buf(buf + j, t3.z);
    }
    
};

typedef ReallocableArray<triangle *> triangle_l;

class space_t : public space {
public:
    triangle_l my_triangles;
    
    space_t()
    { }
    
    virtual ~space_t()
    { erase_all(); }
    
    void add_new_triangle(triangle *t)
    {
        my_triangles.Append(t);
    }
    
    virtual void print_datas()
    {
        fprintf(stderr, "Space of triangle.\n");
                for (int i=0; i< my_triangles.NumElements(); i++) {
                  triangle* cur = my_triangles[i];
          cur->print_datas();
                }
    }
    
    virtual void erase_all()
    { my_triangles.DeleteAll(); }
    
    virtual void send_datas(Socket *s)
    {
        unsigned short size;
        int j;
        char *buf;
        us_t l;
        
        size = 73 * my_triangles.NumElements();
        buf = (char *) malloc(sizeof(char) * (3 + size));
        if(buf == NULL)
            return;
        l.s = htons(size);
        buf[0] = TYPE_TRIANGLE;
        for(int k=0; k<2; k++)
            buf[k+1] = l.c[k];
        j=0;
                for (int i=0; i< my_triangles.NumElements(); i++) {
                  triangle* cur = my_triangles[i];
          (cur)->data_to_send(buf + 3 + (73*j));
          j++;
        }
        s->Send(buf, 3 + size);
        free(buf);
    }
    
    virtual void use_new_datas(char *buf, int size)
    {
        int i,nb;
        triangle *t;
        
        nb = size / 73;
        
        erase_all();
        
        for(i=0;i<nb; i++)
        {
            t = new triangle();
            t->receive_data(buf);
            add_new_triangle(t);
            buf+=73;
        }
    }
    
    virtual Point3D *get_force(Point3D *pos, Point3D *res)
    {
        Point3D tmp;
        
        if(res == NULL)
            res = new Point3D;
        
        *res = Point3D(0,0,0);
        
                for (int i=0; i< my_triangles.NumElements(); i++) {
                  triangle* cur = my_triangles[i];
          cur->get_force(pos, &tmp);
          *res += tmp;
        }
        
        return res;
    }
    
    virtual int istouching()
    {
                for (int i=0; i< my_triangles.NumElements(); i++) {
                  triangle* cur = my_triangles[i];
          if(cur->force)
              return 1;
        }
        return 0;
    }
    
    virtual long gethash()
    {
        long res = 0;
                for (int i=0; i< my_triangles.NumElements(); i++) {
                  triangle* cur = my_triangles[i];
          res ^= cur->gethash();
        }
        
        return res;
    }
    
};

#endif

---