Modelling Natural Shapes: (Easter) Eggs 2020

One year ago, I wrote an article about the modelling of the egg shapes, promising at one point to come back on the topics. A next step in studying eggs shapes is to look to real one or a copy of it. A happy occasion for experimenting with the model using three-dimensional graphics and 3d Printing! That is a natural indeed step: take half of the symmetric curve representing the egg shape

$y=T(1+x)^{\frac{\lambda}{1+\lambda}}(1-x)^{\frac{1}{1+\lambda}}$ ,

where $T$ and $\lambda$ are two parameters, and rotate it around the central axis

$\begin{aligned} x'&=&x\\ y' &=&y*cos(\theta) \\ z' &=& y*sin(\theta) \end{aligned}$

Then using the OpenGL and GLUT graphical libraries, it is quite easy to come out with a simple program that does the job. I have added some screenshots of the program that youc an find in the appendix of the article. In the paper by Stoddard et al. [3], the asymmetry ( $A=\lambda-1$ ) and ellipticity ( $E=\frac{1}{T}-1$ ) parameters are used instead to define the charateristic of the analyzed eggs. For $A=0$ and $E=0$ , the egg has a spherical shape, as in the figure below.

If the value of the parameter $E$ increase the egg shape become more elliptic. For example, with A=0 and E=1, we get the egg in the Figure 2.

If the value of the parameter $A$ increase the egg shape become more elliptic. For example, with A=1 and E=0, we get the egg in the Figure 3.

Please refer to my previous article to have information about the connection of the above egg shapes example with those of real bird eggs.

If you liked this article do not forget to press the like button and share it!

BUONA PASQUA – HAPPY EASTER – FROHE OSTERN

APPENDIX

This is the source code in C++ of the program. It use OpenGl library and the GLUT library for the GUI. The code can be compiled under MacOS using the comand:

g++ eggshapes.cpp -framework opengl -framework cocoa -framework glut -Wno-deprecated -lglui

It should be easy to compile it under Linux, Android or Windows system as well.

/*
 * PROGRAM NAME: Eggshapes 
 * VERSION     : 1.0 
 * DESCRIPTION : 
 *
 * THE BIRDS EGG SHAPE GENERATOR
 *
 * The program implement the Baker equation in 3D 
 * to generate shapes of eggs of different birds 
 *
 *
 * (c) Danilo Roccatano 2019-20
 * 
 */

#include <fstream>
#include <iostream>
#include <iomanip>
#include <math.h>
#include <GLUT/glut.h>
#include <OpenGL/gl.h>
#include <OpenGL/glu.h>
#include <GL/glui.h>

#define DEBUG  -1

using namespace std;

/** These are the live variables passed into GLUI ***/
float  pA=0.5;        // initial value variable T 
float  pE=0.5;        // initial value variable lambda 
int   np=100  ;      // Number of points for the egg shape 
int   nrot=100  ;    // Number of rotation along the egg axis 
int   main_window;

GLint spin= 0;
GLint spin1=0;
GLint egg= 0;
GLUI_Spinner    *spi[11];

//Viewer options (GluLookAt)

float fovy = 60.0, aspect = 1.0, zNear = 1.0, zFar = 100.0;

//Mouse modifiers

float depth = -1;
float scale=2.5;
float psi=0, theta=0;
float downX, downY;
bool leftButton = false, middleButton = false, rightButton=false;
int  showAxis;
int  showLight=1,showLight1=1;
int  sqstripes=-1;
int  chir;

int DrawEgg();
void mouse ();
void atome(GLfloat []);
void drawAxis();
void mouseCallback(int button, int state, int x, int y);
void motionCallback(int x, int y);


/***************************************** myGlutIdle() ***********/

void myGlutIdle( void )
{
    /* According to the GLUT specification, the current window is 
       undefined during an idle callback.  So we need to explicitly change
       it if necessary */
    if ( glutGetWindow() != main_window )
        glutSetWindow(main_window);

    glutPostRedisplay();
}

// ============ Some 3D math routines

//------ Returns the cross product of 2 vectors

void CrossProduct (double M[], double N[], double CS[]) {
    CS[0]=M[1] * N[2] - M[2] * N[1];
    CS[1]=M[2] * N[0] - M[0] * N[2];
    CS[2]=M[0] * N[1] - M[1] * N[0];
}

//------ Returns the length of a vector

double GetVectorLength (double M[]) {
    return sqrt( M[0] * M[0] + M[1] * M[1] + M[2] * M[2] );
}

//------ Returns the sum of two vectors

void   AddVectors (double M[], double N[], double R[]){
    R[0]= M[0] + N[0];
    R[1]= M[1] + N[1];
    R[2]= M[2] + N[2];
}

//------ Returns the vector scaled by the last parameter

void ScaleVector (double M[], double a, double N[]) {
    N[0]= M[0] * a;
    N[1]= M[1] * a; 
    N[2]= M[2] * a;
}

//------ Returns a normalized vector (length = 1)

void NormalizeVector (double M[],double R[3]) {
    double norm = GetVectorLength(M);
    if (norm == 0) norm  =1.0;
    ScaleVector( M, 1./ norm, R );
}

//------ Returns the unit normal vector of a triangle specified by the three
// points P1, P2, and P3.

void FindUnitNormal (double P1[3],double P2[3],double P3[3],double N[3]){
    double D1[3],D2[3];
    double R[3];
    D1[0]=P1[0]-P2[0];
    D1[1]=P1[1]-P2[1];
    D1[2]=P1[2]-P2[2];
    D2[0]=P2[0]-P3[0];
    D2[1]=P2[1]-P3[1];
    D2[2]=P2[2]-P3[2];
    CrossProduct(D1,D2,R);
    NormalizeVector(R,N);
}

//------ Initialization routine

void init () {

    GLfloat  mat_ambient[]  = { 0.25,     0.20725,  0.20725,  0 };
    GLfloat  mat_diffuse[]  = { 1.0,        0.829,    0.829,    0 };
    GLfloat  mat_specular[] = { 0.296648, 0.296648, 0.296648, 0 };
    GLfloat  light_diffuse[]  = { 0.8, 0.8, 0.8, 1.0 };
    GLfloat  light_ambient[]  = { 0.8, 0.8, 0.8, 1.0 };
    GLfloat  light_specular[] = { 0.5, 0.5, 0.5, 1.0 };

    GLfloat  light0_position[] = { -20.0, 20.0, 0.0, 0 };
    GLfloat  light1_position[] = { 10.0, 50.0, 0.0, 0 };

    GLfloat  shininess=0.088 * 128;

    glClearColor( 0  , 0  , 0  , 0   );    // Black background
    glShadeModel(GL_SMOOTH);       // Smooth shading
    glEnable(GL_MULTISAMPLE);      // Enable multisample antialiasing
    glEnable(GL_DEPTH_TEST);       // Enable hidden surface removal
    glEnable(GL_LIGHTING);

    // set light 0

    glEnable(GL_LIGHT0);
    glLightfv( GL_LIGHT0, GL_POSITION, light0_position );
    glLightfv( GL_LIGHT0, GL_DIFFUSE,  light_diffuse );
    glLightfv( GL_LIGHT0, GL_AMBIENT,  light_ambient );
    glLightfv( GL_LIGHT0, GL_SPECULAR, light_specular );

    // set light 1

    glEnable(GL_LIGHT1);
    glLightfv( GL_LIGHT1, GL_POSITION, light1_position );
    glLightfv( GL_LIGHT1, GL_DIFFUSE,  light_diffuse );
    glLightfv( GL_LIGHT1, GL_AMBIENT,  light_ambient );
    glLightfv( GL_LIGHT1, GL_SPECULAR, light_specular );

    // set material

    glMaterialfv( GL_FRONT, GL_AMBIENT,  mat_ambient );
    glMaterialfv( GL_FRONT, GL_DIFFUSE,  mat_diffuse );
    glMaterialfv( GL_FRONT, GL_SPECULAR, light_specular );
    glMaterialf( GL_FRONT, GL_SHININESS, shininess );

    egg=DrawEgg();
}

//------ Draw the scene

void display () {
    glClear( GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT );
    glLoadIdentity();
    gluLookAt( 2, 4, 10, 0, 0, 0, 0, 1, 0 );
    glPushMatrix();
    glScalef( scale, scale, scale );
    //Motion Options

    glTranslatef(0.0, 0.0, -depth);
    glRotatef(-theta, 1.0, 0.0, 0.0);
    glRotatef(psi, 0.0, 1.0, 0.0);


    if (showLight) {
        //        glEnable(GL_LIGHTING);
        glEnable(GL_LIGHT0);
    }
    else
    {
        //        glDisable(GL_LIGHTING);
        glDisable(GL_LIGHT0);
    }

    if (showLight1) {
        //        glEnable(GL_LIGHTING);
        glEnable(GL_LIGHT1);
    }
    else
    {
        //        glDisable(GL_LIGHTING);
        glDisable(GL_LIGHT1);
    }

    glCallList(egg);            // draw the egg 

    glPopMatrix();
    glutSwapBuffers();
}

//------ GLUT Callback called when the window is resized

void reshape (int w, int h)
{
    if (h == 0 || w == 0) return;  //Nothing is visible then, so return
    //Set a new projection matrix

    glViewport( 0, 0, w, h );
    glMatrixMode(GL_PROJECTION);
    glLoadIdentity();
    gluPerspective( 60, h ? w / h : 0, 1, 20 );
    glMatrixMode(GL_MODELVIEW);
    glLoadIdentity();
}

//------ Routine for rotating the scene

void spinDisplay (){

    int TimeNow = glutGet(GLUT_ELAPSED_TIME);
    int WaitUntil=0;

    if ( TimeNow >= WaitUntil ) {
        spin += 1;

        if ( spin > 360 )spin = spin - 360; 
        glutPostRedisplay();
        WaitUntil = TimeNow + 1000 / 25;    // 25 frames/s
    }
}

int DrawEgg(){

    /*
     *       Draw a eggs in 3D using Baker equation 
     */

    // ------ Parameters of the egg 

    double pi=3.141592;

    int in=0;
    double C[np*2][3],C0[np*2][3];
    double NN[3];
    GLuint  egg;


    double dthe= 2*pi/nrot;
    double delta =2.0 / np;
    double lam= pA+1.0;
    double T=1.0/(pE+1.0);
    double ex1=1.0/(1.0+lam);
    double ex2=lam/(1.0+lam);
    double x,x1,x2,y,t;

    //------ Draw the egg 
    egg   = glGenLists(1);
    glNewList( egg, GL_COMPILE );

    for ( double i=0 ; i <= np ; i ++ ) {
        x=-1.0+i * delta;
        x1=1.0+x;
        x2=1.0-x;
        y=T*pow(x1,ex1)*pow(x2,ex2);

        cout << setprecision(6) << fixed;


        // OpenGL display list for the egg   
        for ( int j=0;j<=nrot ;j++ )
        {                      

            t=dthe*j;
            C[j][0] =x;
            C[j][1] =y*cos(t);
            C[j][2] =y*sin(t);
        }

        glBegin(GL_TRIANGLE_STRIP);  // draws the band between the circles C0 and C

        if(sqstripes) {
            double r1= rand() / double(RAND_MAX);
            double r2= rand() / double(RAND_MAX);
            double r3= rand() / double(RAND_MAX);
            GLfloat  mat_ambient[]  = { r1 ,     r2,  r3,  0 };
            glMaterialfv( GL_FRONT, GL_AMBIENT,  mat_ambient );
        } 

        if (i>0) {
            for ( int j=0;j<=nrot ;j++ )
            {
                in=j;

                if(!sqstripes) {
                    double r1= rand() / double(RAND_MAX);
                    double r2= rand() / double(RAND_MAX);
                    double r3= rand() / double(RAND_MAX);
                    GLfloat  mat_ambient[]  = { r1 ,     r2,  r3,  0 };
                    glMaterialfv( GL_FRONT, GL_AMBIENT,  mat_ambient );
                }

                FindUnitNormal (C0[in],C0[j+1],C[in],NN);
                glNormal3d(NN[0],NN[1],NN[2]);

                glVertex3f( C0[j][0],C0[j][1],C0[j][2]  );
                glVertex3f( C[j][0],C[j][1],C[j][2] );
                glVertex3f( C0[j+1][0],C0[j+1][1],C0[j+1][2]);
                glVertex3f( C[j+1][0],C[j+1][1],C[j+1][2]);

            }
        }
        for ( int j=0;j<=nrot ;j++ )
        {
            C0[j][0]=C[j][0];
            C0[j][1]=C[j][1];
            C0[j][2]=C[j][2];
        }

        glEnd();
    }
    //Connect the last with the first 

    glEndList();

    return egg;
}

//------ GLUT callback for the mouse

void  mouse (GLint button, GLint state, GLdouble x, GLdouble y )
{
    if ( button == GLUT_LEFT_BUTTON ) {
        if ( state == GLUT_DOWN ) glutIdleFunc( spinDisplay )   ;
    }
    else if ( button == GLUT_RIGHT_BUTTON ) {
        //        if ( state == GLUT_DOWN ) glutIdleFunc(undef) ;
    }
}

/* Callbacks */
void mouseCallback(int button, int state, int x, int y)
{
    downX = x; downY = y;
    leftButton = ((button == GLUT_LEFT_BUTTON) && (state == GLUT_DOWN));
    rightButton = ((button == GLUT_RIGHT_BUTTON) && (state == GLUT_DOWN));
    middleButton = ((button == GLUT_MIDDLE_BUTTON) &&  (state == GLUT_DOWN));
}

void motionCallback(int x, int y)
{
    if (leftButton) //Rotate
    {
        psi += (x-downX)/4.0;
        theta += (downY-y)/4.0;
    }
    if (rightButton) //Scale
    {
        if (depth + (downY - y)/10.0 < zFar-10 && depth + (downY - y)/10.0 > zNear+3)
            depth += (downY - y)/10.0;
    }
    downX = x;
    downY = y;

    glutPostRedisplay();
}
void control_cb( int control )
{
    egg=DrawEgg();
}


int main (int argc, char **argv)
{
    //Initialize GLUT
    glutInit(&argc, argv);
    //double buffering used to avoid flickering problem in animation
    glutInitDisplayMode(
            GLUT_DOUBLE         // Double buffering
            | GLUT_RGB            // RGB color mode
            | GLUT_DEPTH          // Hidden surface removal
            | GLUT_MULTISAMPLE    // Multisample antialiasing
            );

    // window size
    glutInitWindowSize( 600, 600 );
    // create the window 
    main_window = glutCreateWindow("Egg shapes");
    init();
    //Assign  the function used in events
    glutDisplayFunc(display);

    glutReshapeFunc( reshape );
    //    glutMouseFunc( mouse );
    // motion
    glutMouseFunc(mouseCallback);
    glutMotionFunc(motionCallback);

    /*         Here's the GLUI code         */
    /****************************************/

    GLUI *glui = GLUI_Master.create_glui( "Command & Control",300,500 );
    glui->add_statictext("Parameters for aperture shape");
    spi[0]= new GLUI_Spinner( glui, "Asymmetric parameter (A):", &pA,200,control_cb);
    spi[0] ->set_float_limits( -10, 10. );
    spi[0]->set_speed(0.05);
    spi[1]= new GLUI_Spinner( glui, "Ellipticity parameter (E):", &pE,200,control_cb );
    spi[1]->set_float_limits( -10, 10. );
    spi[1]->set_speed(0.05);

    glui->add_separator();
    new GLUI_Checkbox( glui, "Show Lights 0     ", &showLight);
    new GLUI_Checkbox( glui, "Show Lights 1    ", &showLight1);

    glui->add_separator();
    new GLUI_Checkbox( glui, "Square/Stripes     ", &sqstripes,0,control_cb);

    glui->add_separator();  
    (new GLUI_Spinner( glui, "Scaling:", &scale,206, control_cb )) ->set_float_limits( 0, 20. );

    glui->add_separator();

    glui->add_button( "Quit", 0,(GLUI_Update_CB)exit );
    glui->set_main_gfx_window( main_window );

    /* We register the idle callback with GLUI, *not* with GLUT */
    GLUI_Master.set_glutIdleFunc( myGlutIdle );    
    //Let start glut loop
    glutMainLoop();
    return 0;
}