source: cpp/frams/genetics/fF/conv_fF.cpp @ 740

Last change on this file since 740 was 740, checked in by Maciej Komosinski, 6 years ago

Added support for developmental phases ("checkpoints") to converters f9->f0 and fF->f0

  • Property svn:eol-style set to native
File size: 8.0 KB
RevLine 
[286]1// This file is a part of Framsticks SDK.  http://www.framsticks.com/
[736]2// Copyright (C) 1999-2018  Maciej Komosinski and Szymon Ulatowski.
[286]3// See LICENSE.txt for details.
[140]4
5#include "conv_fF.h"
6#include "fF_genotype.h"
7#include <common/nonstd_stl.h>
[185]8#include <common/Convert.h>
[140]9
[177]10GenoConv_fF0::GenoConv_fF0()
[176]11{
[667]12        name = "10-parameter Foraminifera encoding";
[176]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];
[178]18        precompute_cos_and_sin();
[140]19}
20
[177]21GenoConv_fF0::~GenoConv_fF0()
[176]22{
23        delete[] cosines;
24        delete[] sines;
[174]25}
[140]26
[256]27Part *GenoConv_fF0::addNewPart(Model *m, const fF_chamber3d* c)
28{
29        Part *part = m->addNewPart(Part::SHAPE_ELLIPSOID);
30        part->p = Pt3D(c->centerX, c->centerY, c->centerZ);
31        Pt3D hole = Pt3D(c->holeX, c->holeY, c->holeZ);
32        Orient o;
33        o.lookAt(part->p - hole);
34        part->setOrient(o);
35        return part;
36}
37
[736]38SString GenoConv_fF0::convert(SString &in, MultiMap *map, bool using_checkpoints)
[176]39{
40        fF_growth_params gp;
[348]41        if (!gp.load(in.c_str())) //invalid input genotype?
[176]42                return ""; //so we return an invalid f0 genotype
[140]43
[176]44        Model m;
[740]45        m.open(using_checkpoints);
[177]46
47        m.vis_style = "foram"; //dedicated visual look for Foraminifera
48
[176]49        fF_chamber3d **chambers = new fF_chamber3d*[gp.number_of_chambers];
[177]50        for (int i = 0; i < gp.number_of_chambers; i++)
[667]51                createSphere(i, chambers, gp.radius0x, gp.radius0y, gp.radius0z, gp.translation, gp.angle1, gp.angle2, gp.scalex, gp.scaley, gp.scalez);
[174]52
[256]53        Part *p1 = addNewPart(&m, chambers[0]);
[667]54        p1->scale = Pt3D(gp.radius0x, gp.radius0y, gp.radius0z); //size of the initial chamber
[740]55        m.checkpoint();
[256]56        for (int i = 1; i < gp.number_of_chambers; i++)
[177]57        {
[256]58                Part *p2 = addNewPart(&m, chambers[i]);
[176]59                p2->scale = p1->scale.entrywiseProduct(Pt3D(gp.scalex, gp.scaley, gp.scalez)); //each part's scale is its predecessor's scale * scaling
[544]60                m.addNewJoint(p1, p2, Joint::SHAPE_FIXED); //all parts must be connected
[740]61                m.checkpoint();
[256]62                p1 = p2;
[176]63        }
[174]64
[176]65        for (int i = 0; i < gp.number_of_chambers; i++)
66                delete chambers[i];
[256]67        delete[] chambers;
[174]68
[176]69        m.close();
[544]70        return m.getF0Geno().getGenes();
[140]71}
[174]72
[667]73void GenoConv_fF0::createSphere(int which, fF_chamber3d **chambers, double radius0x, double radius0y, double radius0z, double translation, double alpha_, double gamma_, double kx_, double ky_, double kz_)
[176]74{
[667]75        chambers[which] = new fF_chamber3d(0.0, 0.0, 0.0,
76                radius0x, radius0y, radius0z, radius0x * kx_, 0.0, 0.0,
77                radius0x * translation, 0.0, 0.0, 0.0, 0.0);
[176]78        if (which == 0)
[667]79                chambers[which]->points = generate_points(chambers[which]);
80        if (which > 0)
81        {
82                chambers[which]->radius_x = get_radius(chambers[which - 1]->radius_x, kx_, chambers[0]->radius_x);
83                chambers[which]->radius_y = get_radius(chambers[which - 1]->radius_y, ky_, chambers[0]->radius_y);
84                chambers[which]->radius_z = get_radius(chambers[which - 1]->radius_z, kz_, chambers[0]->radius_z);
85
[176]86                /* new growth vector length */
[667]87                double len = chambers[which]->radius_y * translation;
88                double max_radius = fF_TOO_MUCH * chambers[which]->radius_y;
89                if (fabs(len) > (max_radius))
90                        len = ((len < 0) ? (-1) : 1) * max_radius;
91                if (len == 0)
[176]92                        len = -0.0000001;
[174]93
[176]94                /* aperture of the previous chamber */
95                double pzx = chambers[which - 1]->holeX;
96                double pzy = chambers[which - 1]->holeY;
97                double pzz = chambers[which - 1]->holeZ;
[174]98
[176]99                //center of the previous chamber
100                double pcx = chambers[which - 1]->centerX;
101                double pcy = chambers[which - 1]->centerY;
[319]102                //double pcz = chambers[which - 1]->centerZ; //not used
[174]103
[176]104                /* aperture of the next to last chamber */
105                double ppx;
106                double ppy;
[319]107                //double ppz; //not used
[174]108
[177]109                if (which == 1)
110                {
[176]111                        ppx = pcx;
112                        ppy = pcy;
[319]113                        //ppz = pcz;
[176]114                }
[177]115                else
116                {
[176]117                        ppx = chambers[which - 2]->holeX;
118                        ppy = chambers[which - 2]->holeY;
[319]119                        //ppz = chambers[which - 2]->holeZ;
[176]120                }
[174]121
[176]122                double pzxprim = pzx - ppx;
123                double pzyprim = pzy - ppy;
124                double angle;
[174]125
[185]126                angle = Convert::atan_2(pzyprim, pzxprim);
[176]127                double alpha = angle - alpha_;
[174]128
[176]129                double gamma = chambers[which - 1]->phi + gamma_;
[174]130
[176]131                double wx = len * cos(alpha);
132                double wy = len * sin(alpha);
133                double wz = len * sin(alpha) * sin(gamma);
[174]134
[176]135                /*center of the new sphere*/
136                double x = pzx + wx;
137                double y = pzy + wy;
138                double z = pzz + wz;
[174]139
[667]140                chambers[which]->centerX = x;
141                chambers[which]->centerY = y;
142                chambers[which]->centerZ = z;
143                chambers[which]->vectorTfX = wx;
144                chambers[which]->vectorTfY = wy;
145                chambers[which]->vectorTfZ = wz;
146                chambers[which]->beta = alpha;
147                chambers[which]->phi = gamma;
[174]148
[667]149                chambers[which]->points = generate_points(chambers[which]);
150                search_hid(which, chambers);
[176]151                int pun;
[667]152                pun = find_hole(which, pzx, pzy, pzz, chambers);
153                if (pun < 0) //should never happen
154                {
155                        logPrintf("GenoConv_fF0", "createSphere", LOG_ERROR, "find_hole(%d) returned %d", which, pun);
156                        pun = 0;
157                }
158                chambers[which]->holeX = chambers[which]->points[pun].x;
159                chambers[which]->holeY = chambers[which]->points[pun].y;
160                chambers[which]->holeZ = chambers[which]->points[pun].z;
161        }
162}
[174]163
[667]164double GenoConv_fF0::get_radius(double prev_radius, double scale, double start_radius)
165{
166        double radius = prev_radius * scale;
167        double min_radius = fF_TOO_LITTLE*start_radius;
168        if (radius < min_radius) {
169                radius = min_radius;
[176]170        }
[667]171
172        return radius;
[174]173}
174
[178]175void GenoConv_fF0::precompute_cos_and_sin()
[176]176{
[667]177        double angle = M_PI * 2 / fF_LATITUDE_NUM;
178        for (int i = 0; i < fF_LATITUDE_NUM; i++)
[176]179        {
[667]180                cosines[i] = cos(i * angle);
181                sines[i] = sin(i * angle);
[176]182        }
[174]183}
184
[667]185fF_point* GenoConv_fF0::generate_points(fF_chamber3d *chamber)
[176]186{
[667]187        double cenx = chamber->centerX;
188        double ceny = chamber->centerY;
189        double cenz = chamber->centerZ;
[174]190
[667]191        double rx = chamber->radius_x;
192        double ry = chamber->radius_y;
193        double rz = chamber->radius_z;
[174]194
[178]195        fF_point *points = new fF_point[fF_SIZE];
[174]196
[177]197        for (int i = 0; i < fF_LONGITUDE_NUM; i++)
198        {
199                double y = ceny + ry * cosines[i];
[174]200
[177]201                for (int j = 0; j < fF_LATITUDE_NUM; j++)
202                {
203                        double x = cenx + rx * cosines[j] * sines[i];
204                        double z = cenz + rz * sines[j] * sines[i];
[178]205                        fF_point &p = points[(i * fF_LATITUDE_NUM) + j];
206                        p.x = x;
207                        p.y = y;
208                        p.z = z;
209                        p.inside = false;
[177]210                }
211        }
[176]212        return points;
213
[174]214}
215
[667]216template<typename T> T Square(T x) { return x * x; }
217
[176]218double GenoConv_fF0::dist(double x1, double y1, double z1, double x2, double y2, double z2)
219{
[667]220        return sqrt(Square(x2 - x1) + Square(y2 - y1) + Square(z2 - z1));
[174]221}
222
[667]223void GenoConv_fF0::search_hid(int nr, fF_chamber3d **chambers)
[176]224{
[177]225        for (int i = 0; i < nr; i++)
226        {
[667]227                fF_chamber3d *chamber = chambers[i];
[174]228
[667]229                double rx_sq = Square(chamber->radius_x);
230                double ry_sq = Square(chamber->radius_y);
231                double rz_sq = Square(chamber->radius_z);
[174]232
[177]233                for (int j = 0; j < fF_AMOUNT; j++)
234                {
[667]235                        fF_point &p = chambers[nr]->points[j];
[174]236
[667]237                        double upx = Square(p.x - chamber->centerX);
238                        double upy = Square(p.y - chamber->centerY);
239                        double upz = Square(p.z - chamber->centerZ);
[174]240
[667]241                        double expx = upx / rx_sq;
242                        double expy = upy / ry_sq;
243                        double expz = upz / rz_sq;
[174]244
[667]245                        double result = expx + expy + expz;
[174]246
[177]247                        if (result < fF_THICK_RATIO)
248                        {
[178]249                                p.inside = true;
[176]250                        }
251                }
252        }
[174]253}
254
[667]255int GenoConv_fF0::find_hole(int which, double x, double y, double z, fF_chamber3d **chambers)
[176]256{
[177]257        int found = -1;
258        double distsq_found;
[174]259
[177]260        for (int i = 0; i < fF_AMOUNT; i++)
261        {
[178]262                fF_point &p = chambers[which]->points[i];
263                if (!p.inside) //it is not inside another chamber
[176]264                {
[667]265                        double distancesq = Square(p.x - x) + Square(p.y - y) + Square(p.z - z);
[177]266                        if (found < 0)
267                        {
[176]268                                found = i;
[177]269                                distsq_found = distancesq;
[176]270                        }
[177]271                        if (distancesq < distsq_found)
272                        {
273                                if (which != 0)
274                                {
[176]275                                        bool good = true;
[177]276                                        for (int j = 0; j < which && good; j++)
277                                        {
[667]278                                                fF_chamber3d *chamber = chambers[j];
[174]279
[667]280                                                double rx_sq = Square(chamber->radius_x);
281                                                double ry_sq = Square(chamber->radius_y);
282                                                double rz_sq = Square(chamber->radius_z);
[174]283
[667]284                                                double upx = Square(p.x - chamber->centerX);
285                                                double upy = Square(p.y - chamber->centerY);
286                                                double upz = Square(p.z - chamber->centerZ);
[174]287
[667]288                                                double expx = upx / rx_sq;
289                                                double expy = upy / ry_sq;
290                                                double expz = upz / rz_sq;
[174]291
[667]292                                                double result = expx + expy + expz;
[177]293                                                if (result < 1.0)
294                                                {
295                                                        good = false;
[176]296                                                }
297                                        }
[177]298                                        if (good)
299                                        {
[176]300                                                found = i;
[177]301                                                distsq_found = distancesq;
[176]302                                        }
303                                }
304                        }
305                }
306        }
[174]307
[177]308        return found;
[740]309}
Note: See TracBrowser for help on using the repository browser.