Ignore:
Timestamp:
03/15/14 02:30:20 (11 years ago)
Author:
Maciej Komosinski
Message:

Renamed #defines that conflicted with system-wide defines for Windows (added fF_ prefix), updated ranges for scaling parameters, updated the simplest genotype, formatted sources

File:
1 edited

Legend:

Unmodified
Added
Removed

  • cpp/frams/genetics/fF/conv_fF.cpp

    r174 r176  
    88#include <common/nonstd_stl.h>
    99
    10 GenoConv_fF0::GenoConv_fF0() {
    11     name = "7-value Foraminifera encoding";
    12     in_format = 'F';
    13     out_format = '0';
    14     mapsupport = 0;
    15     cosines = new double[LATITUDE_NUM];
    16     sines = new double[LATITUDE_NUM]; 
    17     fill_cos_and_sin();
    18 }
    19 
    20 GenoConv_fF0::~GenoConv_fF0() {
    21     delete []cosines;
    22     delete []sines;
    23 }
    24 
    25 SString GenoConv_fF0::convert(SString &in, MultiMap *map) {
    26     fF_growth_params gp;
    27     if (!gp.load(in)) //invalid input genotype?
    28         return ""; //so we return an invalid f0 genotype
    29 
    30     double div_radius_length = 1;//div_radius_lenght=1 or kx=ky=kz=1
    31     double radius = 1;
    32 
    33     Model m;
    34     m.open();
    35     // subsequent parts (chambers) are placed relative to the previous part's orientation and location
    36     Part *p1, *p2;
    37 
    38     fF_chamber3d **chambers = new fF_chamber3d*[gp.number_of_chambers];
    39 
    40     for (int i = 0; i < gp.number_of_chambers; i++) {
    41         createSphere(i, chambers, radius, div_radius_length, gp.translation, gp.angle1, gp.angle2, gp.scalex, gp.scaley, gp.scalez);
    42     }
    43 
    44     p1 = m.addNewPart(Part::SHAPE_ELLIPSOID);
    45     p1->p = Pt3D(chambers[0]->centerX, chambers[0]->centerY, chambers[0]->centerZ);
    46 
    47 
    48     for (int i = 1; i < gp.number_of_chambers; i++, p1 = p2) {
    49         p2 = m.addNewPart(Part::SHAPE_ELLIPSOID);
    50         p2->scale = p1->scale.entrywiseProduct(Pt3D(gp.scalex, gp.scaley, gp.scalez)); //each part's scale is its predecessor's scale * scaling
    51 
    52         p2->p = Pt3D(chambers[i]->centerX, chambers[i]->centerY, chambers[i]->centerZ);
    53 
    54         m.addNewJoint(p1, p2, Joint::SHAPE_SOLID); //all parts must be connected
    55     }
    56 
    57     for (int i = 0; i < gp.number_of_chambers; i++)
    58         delete chambers[i];
    59     delete []chambers;
    60 
    61     m.close();
    62     return m.getF0Geno().getGene();
     10GenoConv_fF0::GenoConv_fF0()
     11{
     12        name = "7-value Foraminifera encoding";
     13        in_format = 'F';
     14        out_format = '0';
     15        mapsupport = 0;
     16        cosines = new double[fF_LATITUDE_NUM];
     17        sines = new double[fF_LATITUDE_NUM];
     18        fill_cos_and_sin();
     19}
     20
     21GenoConv_fF0::~GenoConv_fF0()
     22{
     23        delete[] cosines;
     24        delete[] sines;
     25}
     26
     27SString GenoConv_fF0::convert(SString &in, MultiMap *map)
     28{
     29        fF_growth_params gp;
     30        if (!gp.load(in)) //invalid input genotype?
     31                return ""; //so we return an invalid f0 genotype
     32
     33        double div_radius_length = 1;//div_radius_length=1 or kx=ky=kz=1
     34        double radius = 1;
     35
     36        Model m;
     37        m.open();
     38        // subsequent parts (chambers) are placed relative to the previous part's orientation and location
     39        Part *p1, *p2;
     40
     41        fF_chamber3d **chambers = new fF_chamber3d*[gp.number_of_chambers];
     42
     43        for (int i = 0; i < gp.number_of_chambers; i++) {
     44                createSphere(i, chambers, radius, div_radius_length, gp.translation, gp.angle1, gp.angle2, gp.scalex, gp.scaley, gp.scalez);
     45        }
     46
     47        p1 = m.addNewPart(Part::SHAPE_ELLIPSOID);
     48        p1->p = Pt3D(chambers[0]->centerX, chambers[0]->centerY, chambers[0]->centerZ);
     49
     50
     51        for (int i = 1; i < gp.number_of_chambers; i++, p1 = p2) {
     52                p2 = m.addNewPart(Part::SHAPE_ELLIPSOID);
     53                p2->scale = p1->scale.entrywiseProduct(Pt3D(gp.scalex, gp.scaley, gp.scalez)); //each part's scale is its predecessor's scale * scaling
     54
     55                p2->p = Pt3D(chambers[i]->centerX, chambers[i]->centerY, chambers[i]->centerZ);
     56
     57                m.addNewJoint(p1, p2, Joint::SHAPE_SOLID); //all parts must be connected
     58        }
     59
     60        for (int i = 0; i < gp.number_of_chambers; i++)
     61                delete chambers[i];
     62        delete[]chambers;
     63
     64        m.close();
     65        return m.getF0Geno().getGene();
    6366}
    6467
    6568void GenoConv_fF0::createSphere(int which, fF_chamber3d **chambers, double radius_, double div_radius_length_, double div_vector_length_,
    66         double alpha_, double gamma_, double kx_, double ky_, double kz_)
    67 {   
    68     chambers[which]=new fF_chamber3d(0.0f, 0.0f, 0.0f,
    69             (float) radius_, (float) radius_ * (float) kx_, 0.0f, 0.0f,
    70             (float) (radius_ * div_vector_length_), 0.0f, 0.0f, 0.0f, 0.0f);
    71     if (which == 0)
    72         chambers[which]->points = generate_points(chambers[which], which, kx_, ky_, kz_);
    73     if (which > 0) {
    74         /* old radius */
    75         double radiusOld, radius;
    76         radiusOld = chambers[which - 1]->radius;
    77         radius = div_radius_length_ * radiusOld;
    78         /* new growth vector length */
    79         double len = radius * div_vector_length_;
    80         if (radius < TOO_LITTLE) {
    81             radius = TOO_LITTLE;
    82             if (fabs(len) > (TOO_MUCH * radius)) {
    83                 len = ((len < 0) ? (-1) : 1) * TOO_MUCH * radius;
    84             }
    85         }
    86         if (len == 0) {
    87             len = -0.0000001;
    88         }
    89 
    90         /* aperture of the previous chamber */
    91         double pzx = chambers[which - 1]->holeX;
    92         double pzy = chambers[which - 1]->holeY;
    93         double pzz = chambers[which - 1]->holeZ;
    94 
    95         //center of the previous chamber
    96         double pcx = chambers[which - 1]->centerX;
    97         double pcy = chambers[which - 1]->centerY;
    98         double pcz = chambers[which - 1]->centerZ;
    99 
    100         /* aperture of the next to last chamber */
    101         double ppx;
    102         double ppy;
    103         double ppz;
    104 
    105         if (which == 1) {
    106             ppx = pcx;
    107             ppy = pcy;
    108             ppz = pcz;
    109         } else {
    110             ppx = chambers[which - 2]->holeX;
    111             ppy = chambers[which - 2]->holeY;
    112             ppz = chambers[which - 2]->holeZ;
    113         }
    114 
    115         double pzxprim = pzx - ppx;
    116         double pzyprim = pzy - ppy;
    117         double angle;
    118 
    119         angle = atan2(pzyprim, pzxprim);
    120         double alpha = angle - alpha_;
    121 
    122 
    123         double gamma = chambers[which - 1]->phi + gamma_;
    124 
    125         /* x */
    126         double wx = len * cos(alpha);
    127         /* y */
    128         double wy = len * sin(alpha);
    129         /* y */
    130         double wz = len * sin(alpha) * sin(gamma);
    131 
    132         /*center of the new sphere*/
    133         double x = pzx + wx;
    134         double y = pzy + wy;
    135         double z = pzz + wz;
    136 
    137         chambers[which]->centerX = (float) x;
    138         chambers[which]->centerY = (float) y;
    139         chambers[which]->centerZ = (float) z;
    140         chambers[which]->radius= (float) radius;
    141         chambers[which]->vectorTfX = (float) wx;
    142         chambers[which]->vectorTfY = (float) wy;
    143         chambers[which]->vectorTfZ = (float) wz;
    144         chambers[which]->beta = (float) alpha;
    145         chambers[which]->phi = (float) gamma;
    146 
    147         chambers[which]->points = generate_points(chambers[which], which, kx_,  ky_,  kz_);
    148         search_hid(which, chambers, kx_,  ky_,  kz_);
    149         int pun;
    150         pun = find_hole(which, pzx, pzy, pzz, chambers, kx_,  ky_,  kz_);
    151 
    152         chambers[which]->holeX = (float) chambers[which]->points[pun][0];
    153         chambers[which]->holeY = (float) chambers[which]->points[pun][1];
    154         chambers[which]->holeZ = (float) chambers[which]->points[pun][2];
    155     }
    156 }
    157 
    158 void GenoConv_fF0::fill_cos_and_sin() {
    159     int i;
    160     double pi = acos(-1.0);
    161     double angle = pi / (((double) LATITUDE_NUM)*0.5);
    162     for (i = 0; i < LATITUDE_NUM; i++) {
    163         cosines[i] = cos((double) i * angle);
    164         sines[i] = sin((double) i * angle);
    165     }
    166 }
    167 
    168 double** GenoConv_fF0::generate_points(fF_chamber3d *chamber, int which, double kx_, double ky_, double kz_) {
    169     float radius = chamber->radius;
    170     float cenx = chamber->centerX;
    171     float ceny = chamber->centerY;
    172     float cenz = chamber->centerZ;
    173 
    174     double maxX = 0;
    175     double maxY = 0;
    176     double minX = 0;
    177     double minY = 0;
    178     double minZ = 0;
    179 
    180     double kx = 1;
    181     double ky = 1;
    182     double kz = 1;
    183 
    184     if (which > 0) {
    185         for (int kt = 1; kt < (which + 1); kt++) {
    186             kx = kx * kx_;
    187             ky = ky * ky_;
    188             kz = kz * kz_;
    189         }
    190     }
    191 
    192     int i, j;
    193     double x, y, z;
    194    
    195     double **points = new double*[SIZE];
    196     for (int i = 0; i < SIZE; i++) {
    197         points[i] = new double[4];
    198     }
    199    
    200     for (i = 0; i < LONGITUDE_NUM; i++) {
    201         if (kx_ == 1 && ky_ == 1 && kz_ == 1) {
    202             y = ceny + radius * cosines[i];
    203         }
    204         else {
    205             y = ceny + ky * cosines[i];
    206         }
    207         for (j = 0; j < LATITUDE_NUM; j++) {
    208             if (kx_ == 1 && ky_ == 1 && kz_ == 1) {
    209                 points[(i * LATITUDE_NUM) + j][0] = x = cenx + radius * cosines[j] * sines[i];
    210                 points[(i * LATITUDE_NUM) + j][1] = y;
    211                 points[(i * LATITUDE_NUM) + j][2] = z = cenz + radius * sines[j] * sines[i];
    212             } else {
    213                 points[(i * LATITUDE_NUM) + j][0] = x = cenx + kx * cosines[j] * sines[i];
    214                 points[(i * LATITUDE_NUM) + j][1] = y;
    215                 points[(i * LATITUDE_NUM) + j][2] = z = cenz + kz * sines[j] * sines[i];
    216             }
    217 
    218             points[(i * LATITUDE_NUM) + j][3] = 1.0;
    219             if (x < minX) minX = x;
    220             if (x > maxX) maxX = x;
    221             if (y < minY) minY = y;
    222             if (y > maxY) maxY = y;
    223 
    224             if (z < minZ) minZ = z;
    225         };
    226     };
    227     return points;
    228 
    229 }
    230 
    231 double GenoConv_fF0::dist(double x1, double y1, double z1, double x2, double y2, double z2) {
    232     return sqrt((x2 - x1)*(x2 - x1) + (y2 - y1)*(y2 - y1) + (z2 - z1)*(z2 - z1));
    233 }
    234 
    235 void GenoConv_fF0::search_hid(int nr, fF_chamber3d **spheres, double kx_, double ky_, double kz_) {
    236 
    237     int i, j;
    238     if (nr != 0) {
    239         for (i = 0; i < nr; i++) {
    240             for (j = 0; j < AMOUNT; j++) {
    241                 double X = spheres[nr]->points[j][0];
    242                 double Y = spheres[nr]->points[j][1];
    243                 double Z = spheres[nr]->points[j][2];
    244 
    245                 double srX0 = spheres[i]->centerX;
    246                 double srY0 = spheres[i]->centerY;
    247                 double srZ0 = spheres[i]->centerZ;
    248 
    249                 double a2;
    250                 double b2;
    251                 double c2;
    252 
    253                 if (kx_ != 1) {
    254                     a2 = (kx_ * kx_);
    255                 } else {
    256                     a2 = (spheres[i]->radius * spheres[i]->radius);
    257                 }
    258 
    259                 if (ky_ != 1) {
    260                     b2 = (ky_ * ky_);
    261 
    262                 } else {
    263                     b2 = (spheres[i]->radius * spheres[i]->radius);
    264                 }
    265 
    266                 c2 = (kz_ * spheres[i]->radius) * (kz_ * spheres[i]->radius);
    267 
    268                 double up1 = (X - srX0) * (X - srX0);
    269                 double up2 = (Y - srY0) * (Y - srY0);
    270                 double up3 = (Z - srZ0) * (Z - srZ0);
    271 
    272                 double exp = up1 / a2;
    273                 double exp2 = up2 / b2;
    274                 double exp3 = up3 / c2;
    275 
    276                 double result = exp + exp2 + exp3;
    277 
    278                 if (result < (THICK_RATIO)
    279                         ) {
    280                     spheres[nr]->points[j][3] = 0;
    281                 }
    282             }
    283         }
    284     }
    285 }
    286 
    287 int GenoConv_fF0::find_hole(int which, double x, double y, double z, fF_chamber3d **chambers, double kx_, double ky_, double kz_) {
    288     int i;
    289     double distance;
    290     int found = 0;
    291     double dist_found = 0;
    292     int first = 1;
    293 
    294     for (i = 0; i < AMOUNT; i++) {
    295         if (chambers[which]->points[i][3] != 0) //nie jest wewnatrz inne komory
    296         {
    297             distance = sqrt((chambers[which]->points[i][0] - x)*(chambers[which]->points[i][0] - x) +
    298                     (chambers[which]->points[i][1] - y)*(chambers[which]->points[i][1] - y)+
    299                     (chambers[which]->points[i][2] - z)*(chambers[which]->points[i][2] - z));
    300             if (first != 0) {
    301                 found = i;
    302                 dist_found = distance;
    303                 first = 0;
    304             };
    305             if (distance < dist_found) {
    306                 if (which != 0) {
    307                     bool good = true;
    308                     for (int j = 0; j < which; j++) {
    309                         {
    310                             double X = chambers[which]->points[i][0];
    311                             double Y = chambers[which]->points[i][1];
    312                             double Z = chambers[which]->points[i][2];
    313 
    314                             double srX0 = chambers[j]->centerX;
    315                             double srY0 = chambers[j]->centerY;
    316                             double srZ0 = chambers[j]->centerZ;
    317 
    318                             double a2 = (kx_ * chambers[j]->radius) * (kx_ * chambers[j]->radius);
    319                             double b2 = (ky_ * chambers[j]->radius) * (ky_ * chambers[j]->radius);
    320                             double c2 = (kz_ * chambers[j]->radius) * (kz_ * chambers[j]->radius);
    321 
    322                             double up1 = (X - srX0) * (X - srX0);
    323                             double up2 = (Y - srY0) * (Y - srY0);
    324                             double up3 = (Z - srZ0) * (Z - srZ0);
    325 
    326                             double exp1 = up1 / a2;
    327                             double exp2 = up2 / b2;
    328                             double exp3 = up3 / c2;
    329 
    330                             double result = exp1 + exp2 + exp3;
    331                             if (result < 1.0)
    332                             {
    333                                 good = false;
    334                             }
    335                         }
    336                     }
    337                     if (good) {
    338                         found = i;
    339                         dist_found = distance;
    340                     }
    341                 }
    342             };
    343         };
    344     };
    345 
    346     return (found);
    347 }
     69        double alpha_, double gamma_, double kx_, double ky_, double kz_)
     70{
     71        chambers[which] = new fF_chamber3d(0.0f, 0.0f, 0.0f,
     72                (float)radius_, (float)radius_ * (float)kx_, 0.0f, 0.0f,
     73                (float)(radius_ * div_vector_length_), 0.0f, 0.0f, 0.0f, 0.0f);
     74        if (which == 0)
     75                chambers[which]->points = generate_points(chambers[which], which, kx_, ky_, kz_);
     76        if (which > 0) {
     77                /* old radius */
     78                double radiusOld, radius;
     79                radiusOld = chambers[which - 1]->radius;
     80                radius = div_radius_length_ * radiusOld;
     81                /* new growth vector length */
     82                double len = radius * div_vector_length_;
     83                if (radius < fF_TOO_LITTLE) {
     84                        radius = fF_TOO_LITTLE;
     85                        if (fabs(len) >(fF_TOO_MUCH * radius)) {
     86                                len = ((len < 0) ? (-1) : 1) * fF_TOO_MUCH * radius;
     87                        }
     88                }
     89                if (len == 0) {
     90                        len = -0.0000001;
     91                }
     92
     93                /* aperture of the previous chamber */
     94                double pzx = chambers[which - 1]->holeX;
     95                double pzy = chambers[which - 1]->holeY;
     96                double pzz = chambers[which - 1]->holeZ;
     97
     98                //center of the previous chamber
     99                double pcx = chambers[which - 1]->centerX;
     100                double pcy = chambers[which - 1]->centerY;
     101                double pcz = chambers[which - 1]->centerZ;
     102
     103                /* aperture of the next to last chamber */
     104                double ppx;
     105                double ppy;
     106                double ppz;
     107
     108                if (which == 1) {
     109                        ppx = pcx;
     110                        ppy = pcy;
     111                        ppz = pcz;
     112                }
     113                else {
     114                        ppx = chambers[which - 2]->holeX;
     115                        ppy = chambers[which - 2]->holeY;
     116                        ppz = chambers[which - 2]->holeZ;
     117                }
     118
     119                double pzxprim = pzx - ppx;
     120                double pzyprim = pzy - ppy;
     121                double angle;
     122
     123                angle = atan2(pzyprim, pzxprim);
     124                double alpha = angle - alpha_;
     125
     126
     127                double gamma = chambers[which - 1]->phi + gamma_;
     128
     129                /* x */
     130                double wx = len * cos(alpha);
     131                /* y */
     132                double wy = len * sin(alpha);
     133                /* y */
     134                double wz = len * sin(alpha) * sin(gamma);
     135
     136                /*center of the new sphere*/
     137                double x = pzx + wx;
     138                double y = pzy + wy;
     139                double z = pzz + wz;
     140
     141                chambers[which]->centerX = (float)x;
     142                chambers[which]->centerY = (float)y;
     143                chambers[which]->centerZ = (float)z;
     144                chambers[which]->radius = (float)radius;
     145                chambers[which]->vectorTfX = (float)wx;
     146                chambers[which]->vectorTfY = (float)wy;
     147                chambers[which]->vectorTfZ = (float)wz;
     148                chambers[which]->beta = (float)alpha;
     149                chambers[which]->phi = (float)gamma;
     150
     151                chambers[which]->points = generate_points(chambers[which], which, kx_, ky_, kz_);
     152                search_hid(which, chambers, kx_, ky_, kz_);
     153                int pun;
     154                pun = find_hole(which, pzx, pzy, pzz, chambers, kx_, ky_, kz_);
     155
     156                chambers[which]->holeX = (float)chambers[which]->points[pun][0];
     157                chambers[which]->holeY = (float)chambers[which]->points[pun][1];
     158                chambers[which]->holeZ = (float)chambers[which]->points[pun][2];
     159        }
     160}
     161
     162void GenoConv_fF0::fill_cos_and_sin()
     163{
     164        int i;
     165        double pi = acos(-1.0);
     166        double angle = pi / (((double)fF_LATITUDE_NUM)*0.5);
     167        for (i = 0; i < fF_LATITUDE_NUM; i++)
     168        {
     169                cosines[i] = cos((double)i * angle);
     170                sines[i] = sin((double)i * angle);
     171        }
     172}
     173
     174double** GenoConv_fF0::generate_points(fF_chamber3d *chamber, int which, double kx_, double ky_, double kz_)
     175{
     176        float radius = chamber->radius;
     177        float cenx = chamber->centerX;
     178        float ceny = chamber->centerY;
     179        float cenz = chamber->centerZ;
     180
     181        double maxX = 0;
     182        double maxY = 0;
     183        double minX = 0;
     184        double minY = 0;
     185        double minZ = 0;
     186
     187        double kx = 1;
     188        double ky = 1;
     189        double kz = 1;
     190
     191        if (which > 0) {
     192                for (int kt = 1; kt < (which + 1); kt++) {
     193                        kx = kx * kx_;
     194                        ky = ky * ky_;
     195                        kz = kz * kz_;
     196                }
     197        }
     198
     199        int i, j;
     200        double x, y, z;
     201
     202        double **points = new double*[fF_SIZE];
     203        for (int i = 0; i < fF_SIZE; i++) {
     204                points[i] = new double[4];
     205        }
     206
     207        for (i = 0; i < fF_LONGITUDE_NUM; i++) {
     208                if (kx_ == 1 && ky_ == 1 && kz_ == 1) {
     209                        y = ceny + radius * cosines[i];
     210                }
     211                else {
     212                        y = ceny + ky * cosines[i];
     213                }
     214                for (j = 0; j < fF_LATITUDE_NUM; j++) {
     215                        if (kx_ == 1 && ky_ == 1 && kz_ == 1) {
     216                                points[(i * fF_LATITUDE_NUM) + j][0] = x = cenx + radius * cosines[j] * sines[i];
     217                                points[(i * fF_LATITUDE_NUM) + j][1] = y;
     218                                points[(i * fF_LATITUDE_NUM) + j][2] = z = cenz + radius * sines[j] * sines[i];
     219                        }
     220                        else {
     221                                points[(i * fF_LATITUDE_NUM) + j][0] = x = cenx + kx * cosines[j] * sines[i];
     222                                points[(i * fF_LATITUDE_NUM) + j][1] = y;
     223                                points[(i * fF_LATITUDE_NUM) + j][2] = z = cenz + kz * sines[j] * sines[i];
     224                        }
     225
     226                        points[(i * fF_LATITUDE_NUM) + j][3] = 1.0;
     227                        if (x < minX) minX = x;
     228                        if (x > maxX) maxX = x;
     229                        if (y < minY) minY = y;
     230                        if (y > maxY) maxY = y;
     231
     232                        if (z < minZ) minZ = z;
     233                };
     234        };
     235        return points;
     236
     237}
     238
     239double GenoConv_fF0::dist(double x1, double y1, double z1, double x2, double y2, double z2)
     240{
     241        return sqrt((x2 - x1)*(x2 - x1) + (y2 - y1)*(y2 - y1) + (z2 - z1)*(z2 - z1));
     242}
     243
     244void GenoConv_fF0::search_hid(int nr, fF_chamber3d **spheres, double kx_, double ky_, double kz_)
     245{
     246        int i, j;
     247        if (nr != 0) {
     248                for (i = 0; i < nr; i++) {
     249                        for (j = 0; j < fF_AMOUNT; j++) {
     250                                double X = spheres[nr]->points[j][0];
     251                                double Y = spheres[nr]->points[j][1];
     252                                double Z = spheres[nr]->points[j][2];
     253
     254                                double srX0 = spheres[i]->centerX;
     255                                double srY0 = spheres[i]->centerY;
     256                                double srZ0 = spheres[i]->centerZ;
     257
     258                                double a2;
     259                                double b2;
     260                                double c2;
     261
     262                                if (kx_ != 1) {
     263                                        a2 = (kx_ * kx_);
     264                                }
     265                                else {
     266                                        a2 = (spheres[i]->radius * spheres[i]->radius);
     267                                }
     268
     269                                if (ky_ != 1) {
     270                                        b2 = (ky_ * ky_);
     271
     272                                }
     273                                else {
     274                                        b2 = (spheres[i]->radius * spheres[i]->radius);
     275                                }
     276
     277                                c2 = (kz_ * spheres[i]->radius) * (kz_ * spheres[i]->radius);
     278
     279                                double up1 = (X - srX0) * (X - srX0);
     280                                double up2 = (Y - srY0) * (Y - srY0);
     281                                double up3 = (Z - srZ0) * (Z - srZ0);
     282
     283                                double exp = up1 / a2;
     284                                double exp2 = up2 / b2;
     285                                double exp3 = up3 / c2;
     286
     287                                double result = exp + exp2 + exp3;
     288
     289                                if (result < (fF_THICK_RATIO)
     290                                        ) {
     291                                        spheres[nr]->points[j][3] = 0;
     292                                }
     293                        }
     294                }
     295        }
     296}
     297
     298int GenoConv_fF0::find_hole(int which, double x, double y, double z, fF_chamber3d **chambers, double kx_, double ky_, double kz_)
     299{
     300        int i;
     301        double distance;
     302        int found = 0;
     303        double dist_found = 0;
     304        int first = 1;
     305
     306        for (i = 0; i < fF_AMOUNT; i++) {
     307                if (chambers[which]->points[i][3] != 0) //it is not inside another chamber
     308                {
     309                        distance = sqrt((chambers[which]->points[i][0] - x)*(chambers[which]->points[i][0] - x) +
     310                                (chambers[which]->points[i][1] - y)*(chambers[which]->points[i][1] - y) +
     311                                (chambers[which]->points[i][2] - z)*(chambers[which]->points[i][2] - z));
     312                        if (first != 0) {
     313                                found = i;
     314                                dist_found = distance;
     315                                first = 0;
     316                        }
     317                        if (distance < dist_found) {
     318                                if (which != 0) {
     319                                        bool good = true;
     320                                        for (int j = 0; j < which; j++) {
     321                                                {
     322                                                        double X = chambers[which]->points[i][0];
     323                                                        double Y = chambers[which]->points[i][1];
     324                                                        double Z = chambers[which]->points[i][2];
     325
     326                                                        double srX0 = chambers[j]->centerX;
     327                                                        double srY0 = chambers[j]->centerY;
     328                                                        double srZ0 = chambers[j]->centerZ;
     329
     330                                                        double a2 = (kx_ * chambers[j]->radius) * (kx_ * chambers[j]->radius);
     331                                                        double b2 = (ky_ * chambers[j]->radius) * (ky_ * chambers[j]->radius);
     332                                                        double c2 = (kz_ * chambers[j]->radius) * (kz_ * chambers[j]->radius);
     333
     334                                                        double up1 = (X - srX0) * (X - srX0);
     335                                                        double up2 = (Y - srY0) * (Y - srY0);
     336                                                        double up3 = (Z - srZ0) * (Z - srZ0);
     337
     338                                                        double exp1 = up1 / a2;
     339                                                        double exp2 = up2 / b2;
     340                                                        double exp3 = up3 / c2;
     341
     342                                                        double result = exp1 + exp2 + exp3;
     343                                                        if (result < 1.0)
     344                                                        {
     345                                                                good = false;
     346                                                        }
     347                                                }
     348                                        }
     349                                        if (good) {
     350                                                found = i;
     351                                                dist_found = distance;
     352                                        }
     353                                }
     354                        }
     355                }
     356        }
     357
     358        return (found);
     359}
Note: See TracChangeset for help on using the changeset viewer.