How to show visible part of planar world rendered with 3D perspective on topside 2D minimap?

Prologue

This Q&A is a remake of:

How to highlight 3D perspective camera view (frustrum) on topside 2D minimap?

which was closed (and failed first reopening cycle) due to lack of info and no response of the original author. However I think this is an interesting question though so I decided to Ask&Answer this myself (this time with all the needed specs).

Question

Let assume our world is a uniform rectangular square grid (represented by 2D array of tiles) mapped on a plane (let say plane XY (Z=0.0) for simplicity) and is rendered with perspective projection. like this:

How to map the perspective frustrum (visible part of the map/plane) to the red colored polygonal shape on the minimap ?

To be more universal let assume this as input:

plane (Z=0.0) defined as start point p0 and two basis vectors du,dv which maps the 2D map array of tiles to 3D ...
ModelView matrix and Perspective matrix used

And wanted output:

4 point polygon (on minimap) representing the visible part of out plane

Limitations (to more or less match the original question):

use C++
old style OpenGL (GL,GLU)
no 3th party lib for vector/matrix math

Solution

So what we want is to obtain the 4 intersection points between our plane (Z=0.0) and camera frustrum. So the idea is to cast 4 rays (one for each edge of frustrum) from the camera focal point and simply compute the ray/plane intersection. As the plane is Z=0.0 the intersection point has Z=0.0 too so the intersection is quite easy to compute.

Cast ray for each corner/edge

from camera focal point to screen corner (in screen space)

and convert it to world global coordinates (by reverting perspective and using inverse modelview matrix it is described later). The ray should be in form:
```
p(t) = p + dp*t
```
where p is the focal point and dp is direction vector (does not need to be normalized)
compute the intersection with XY plane (Z=0.0)

As the z=0.0 then:
```
0 = p.z + dp.z*t
t = -p.z/dp.z
```
so we can compute the intersection point directly.
convert 3D intersection points to u,v inside map

for that simple dot product is enough. So if p is our intersection point then:
```
u = dot(p-p0,du)
v = dot(p-p0,dv)
```
where u,v are coordinates in our 2D map array or minimap. In case your u,v are axis aligned then you can use directly (p.x-p0.x,p.y-p0.y) without any dot product

How to convert point p from camera coordinates to global world coordinates:

revert perspective

first obtain perspective matrix parameters
```
double per[16],zNear,zFar,fx,fy;
glGetDoublev(GL_PROJECTION_MATRIX,per);
zFar =0.5*per[14]*(1.0-((per[10]-1.0)/(per[10]+1.0)));
zNear=zFar*(per[10]+1.0)/(per[10]-1.0);
fx=per[0];
fy=per[5];
```
This will give you the frustrums near and far planes and scaling for x,y axises. Now reverting perspective is simply inverting the perspective divide like this:
```
p[1]*=(-p[2]/fy);                   // apply inverse of perspective
p[0]*=(-p[2]/fx);
```
The znear and zfar are needed for casting the rays. For more info see:
- depth buffer got by glReadPixels is always 1
global world coordinates

simply use inverse of ModelView matrix on our p. So first obtain the matrix:
```
double cam[16];
glGetDoublev(GL_MODELVIEW_MATRIX,cam);
```
As inverse you can use my matrix_inv so now the final step is:
```
p = Inverse(cam)*p;
```
but do not forget that p must be homogenuous so (x,y,z,1) for points and (x,y,z,0) for vectors.

Look here if you lack the background knowledge or need vector/matrix math:

Understanding 4x4 homogenous transform matrices

Here Small C++ example of this:

//---------------------------------------------------------------------------
void matrix_mul_vector(double *c,double *a,double *b)
    {
    double q[3];
    q[0]=(a[ 0]*b[0])+(a[ 4]*b[1])+(a[ 8]*b[2])+(a[12]);
    q[1]=(a[ 1]*b[0])+(a[ 5]*b[1])+(a[ 9]*b[2])+(a[13]);
    q[2]=(a[ 2]*b[0])+(a[ 6]*b[1])+(a[10]*b[2])+(a[14]);
    for(int i=0;i<3;i++) c[i]=q[i];
    }
//---------------------------------------------------------------------------
void  matrix_inv(double *a,double *b) // a[16] = Inverse(b[16])
    {
    double x,y,z;
    // transpose of rotation matrix
    a[ 0]=b[ 0];
    a[ 5]=b[ 5];
    a[10]=b[10];
    x=b[1]; a[1]=b[4]; a[4]=x;
    x=b[2]; a[2]=b[8]; a[8]=x;
    x=b[6]; a[6]=b[9]; a[9]=x;
    // copy projection part
    a[ 3]=b[ 3];
    a[ 7]=b[ 7];
    a[11]=b[11];
    a[15]=b[15];
    // convert origin: new_pos = - new_rotation_matrix * old_pos
    x=(a[ 0]*b[12])+(a[ 4]*b[13])+(a[ 8]*b[14]);
    y=(a[ 1]*b[12])+(a[ 5]*b[13])+(a[ 9]*b[14]);
    z=(a[ 2]*b[12])+(a[ 6]*b[13])+(a[10]*b[14]);
    a[12]=-x;
    a[13]=-y;
    a[14]=-z;
    }
//---------------------------------------------------------------------------
void draw_map()
    {
    int i,j;
    double u,v,p[3],dp[3];

    // here 3D view must be already set (modelview,projection)
    glDisable(GL_CULL_FACE);

    // [draw 3D map]
    const int n=30;                 // map size
    double p0[3]={0.0,0.0,0.0};     // map start point
    double du[3]={1.0,0.0,0.0};     // map u step (size of grid = 1.0 )
    double dv[3]={0.0,1.0,0.0};     // map v step (size of grid = 1.0 )
    glColor3f(0.5,0.7,1.0);
    glBegin(GL_LINES);
    for (j=0;j<=n;j++)
        {
        for (i=0;i<3;i++) p[i]=p0[i]+(double(j)*du[i])+(double(0)*dv[i]); glVertex3dv(p);
        for (i=0;i<3;i++) p[i]=p0[i]+(double(j)*du[i])+(double(n)*dv[i]); glVertex3dv(p);
        for (i=0;i<3;i++) p[i]=p0[i]+(double(0)*du[i])+(double(j)*dv[i]); glVertex3dv(p);
        for (i=0;i<3;i++) p[i]=p0[i]+(double(n)*du[i])+(double(j)*dv[i]); glVertex3dv(p);
        }
    glEnd();

    // [compute trapeze points]
    double cam[16],per[16],pt[4][3],zNear,zFar,fx,fy;
    glGetDoublev(GL_PROJECTION_MATRIX,per); // obtain matrices
    glGetDoublev(GL_MODELVIEW_MATRIX,cam);
    matrix_inv(cam,cam);
    zFar =0.5*per[14]*(1.0-((per[10]-1.0)/(per[10]+1.0)));
    zNear=zFar*(per[10]+1.0)/(per[10]-1.0);
    fx=per[0];
    fy=per[5];
    for (j=0;j<4;j++)                       // 4 corners
        {
        for (i=0;i<3;i++) dp[i]=0.0;        // cast ray from camera focus dp
        if (j==0) { p[0]=-1.0; p[1]=-1.0; } // to screen corner p
        if (j==1) { p[0]=-1.0; p[1]=+1.0; }
        if (j==2) { p[0]=+1.0; p[1]=+1.0; }
        if (j==3) { p[0]=+1.0; p[1]=-1.0; }
        p[2]=zNear;                         // start position at screen plane
        p[1]*=(-p[2]/fy);                   // apply inverse of perspective
        p[0]*=(-p[2]/fx);
        // transform to worlds global coordinates
        matrix_mul_vector( p,cam, p);
        matrix_mul_vector(dp,cam,dp);
        // compute intersection of ray and XY plane (z=0) as pt[j] (i exploited the fact that the intersection have z=0.0 for arbitrary plane it would be a bit more complicated)
        for (i=0;i<3;i++) dp[i]=p[i]-dp[i];
        u=p[2]/dp[2];
        if (u<0.0) u=(p[2]-zFar)/dp[2];     // no intersection means "infinite" visibility
        for (i=0;i<3;i++) pt[j][i]=p[i]-(u*dp[i]);
        u=0.0;
        }

    // [draw 2D minimap]
    GLint vp0[4];
    GLint vp1[4]={10,10,150,150};           // minimap position and size ppixels[
    double q0[2]={-1.0,-1.0 };              // minimap start point
    double eu[2]={2.0/double(n),0.0};       // minimap u step
    double ev[2]={0.0,2.0/double(n)};       // minimap v step

    // set 2D view for minimap
    glDisable(GL_DEPTH_TEST);
    glMatrixMode(GL_PROJECTION);
    glPushMatrix();
    glLoadIdentity();
    glMatrixMode(GL_MODELVIEW);
    glPushMatrix();
    glLoadIdentity();
    glGetIntegerv(GL_VIEWPORT,vp0);
    glViewport(vp1[0],vp1[1],vp1[2],vp1[3]);

    glColor3f(0.0,0.0,0.0);                 // clear background
    glBegin(GL_QUADS);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(0)*eu[i])+(double(0)*ev[i]); glVertex2dv(p);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(n)*eu[i])+(double(0)*ev[i]); glVertex2dv(p);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(n)*eu[i])+(double(n)*ev[i]); glVertex2dv(p);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(0)*eu[i])+(double(n)*ev[i]); glVertex2dv(p);
    glEnd();

    glColor3f(0.15,0.15,0.15);              // grid
    glBegin(GL_LINES);
    for (j=0;j<=n;j++)
        {
        for (i=0;i<2;i++) p[i]=q0[i]+(double(j)*eu[i])+(double(0)*ev[i]); glVertex2dv(p);
        for (i=0;i<2;i++) p[i]=q0[i]+(double(j)*eu[i])+(double(n)*ev[i]); glVertex2dv(p);
        for (i=0;i<2;i++) p[i]=q0[i]+(double(0)*eu[i])+(double(j)*ev[i]); glVertex2dv(p);
        for (i=0;i<2;i++) p[i]=q0[i]+(double(n)*eu[i])+(double(j)*ev[i]); glVertex2dv(p);
        }
    glEnd();

    glColor3f(0.5,0.5,0.5);                 // border of minimap
    glLineWidth(2.0);
    glBegin(GL_LINE_LOOP);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(0)*eu[i])+(double(0)*ev[i]); glVertex2dv(p);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(n)*eu[i])+(double(0)*ev[i]); glVertex2dv(p);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(n)*eu[i])+(double(n)*ev[i]); glVertex2dv(p);
    for (i=0;i<2;i++) p[i]=q0[i]+(double(0)*eu[i])+(double(n)*ev[i]); glVertex2dv(p);
    glEnd();
    glLineWidth(1.0);

    // 2D minimap render of the pt[]
    glColor3f(0.7,0.1,0.1);                 // trapeze
    glBegin(GL_LINE_LOOP);
    for (j=0;j<4;j++)
        {
        // get u,v from pt[j]
        for (i=0;i<3;i++) p[i]=pt[j][i]-p0[i];
        for (u=0.0,i=0;i<3;i++) u+=p[i]*du[i];
        for (v=0.0,i=0;i<3;i++) v+=p[i]*dv[i];
        // convert to 2D position and render
        for (i=0;i<2;i++) p[i]=q0[i]+(u*eu[i])+(v*ev[i]); glVertex2dv(p);
        }
    glEnd();

    // restore 3D view
    glMatrixMode(GL_MODELVIEW);
    glPopMatrix();
    glMatrixMode(GL_PROJECTION);
    glPopMatrix();
    glViewport(vp0[0],vp0[1],vp0[2],vp0[3]);
    glEnable(GL_DEPTH_TEST);
    }
//---------------------------------------------------------------------------

And preview:

As you can see we need just matrix*vector multiplication and pseudo inverse matrix functions for this (all others like dot,+,- are really simple and directly encoded as inline code) and both are simple enough to directly implement it in code so no need for GLM or similar lib.

Also I was too lazy to clip the 4 point polygon to minimap size so instead I used glViewport which did it for me.

Here Win32 BDS2006 VCL/C++/OpenGL1.0 Demo:

Demo Source+Binary

Just select slow download and enter the validation code from image. It does not use any 3th party libs other than GL,GLU. The camera is static so just add keyboard/mouse events to your liking. If you want to port this to your environment just mimic the events behavior and ignore the VCL stuff.

The OpenGL init is done based on this:

simple complete GL+VAO/VBO+GLSL+shaders example in C++

I just removed the GLEW,GLSL and VAO stuff from it.