Without delving into your code, I’m going to assume that you build a cube in a three-dimensional space and have some projection transformation that maps this onto a two-dimensional viewport. If so, then you can take advantage of this projection to simplify the problem to one of mapping a rectangle onto the face of your cube. Once you have that mapping, you then apply the same projection to the three-dimensional points of the mapped image, that is, compose the two transformations.
Computing the rectangle-to-cube-face mapping is a fairly simple task. I assume that the image to be mapped is a rectangle of width $W$ and height $H$ with its lower-left corner at the origin and sides parallel to the coordinate axes. Let $V_0(x_0,y_0,z_0)$, $V_1(x_1,y_1,z_1)$, $V_2(x_2,y_2,z_2)$ and $V_3(x_3,y_3,z_3)$ be the vertices of the cube face to which the lower-left, lower-right, upper-left and upper-right corners of the image, respectively, are to be mapped. Note that we must have $V_3=V_1+V_2-V_0$ for these points to define a face of a cube. So, we want an affine transformation that makes the following mappings:$$\begin{align} (0,0)&\mapsto V_0 \\ (W,0)&\mapsto V_1 \\ (0,H)&\mapsto V_2 \\ (W,H)&\mapsto V_1+V_2-V_0. \end{align}$$ The last constraint is redundant, so we won’t be making explicit use of it. A straigforward way to compute this mapping is via homogeneous transformation matrices. Recalling that each column of a transformation matrix is the image of a basis vector, we can construct the matrix of the image-to-cube-face mapping: $$M=\pmatrix{x_0&x_1&x_2\\y_1&y_2&y_3\\z_1&z_2&z_3\\1&1&1}\pmatrix{0&W&0\\0&0&H\\1&1&1}^{-1}=\pmatrix{{x_1-x_0\over W}&{x_2-x_0\over H}&x_0\\{y_1-y_0\over W}&{y_2-y_0\over H}&y_0\\{z_1-z_0\over W}&{z_2-z_0\over H}&z_0\\0&0&1}.$$ The columns of the second matrix in the product on the left are the homogeneous coordinates of the three image corners and the columns of the first matrix are the homogeneous coordinates of the corresponding cube vertices. So, a point $(X,Y)$ on the image is mapped to the point in the cube’s 3-D space with $$x=(x_1-x_0)\frac X W+(x_2-x_0)\frac Y H+x_0$$ (and similarly for $y$ and $z$). You can then apply the view projection to these coordinates to get the resulting coordinates of the image point in the viewpoort.
Note that if you’re applying some transformation to the cube in order to rotate or otherwise move it (or the camera), then you can include that transformation in the pipeline. Compute the image-to-cube-face map once for the cube in its “home” position, and then compose this fixed mapping with the cube transformation to get the mapping of the image onto the transformed cube.
If instead you need to compute the mapping from the image to the viewport directly, the calculations are a bit messier, but still conceptually pretty simple. Here I assume that the projection from the 3-D model space to the viewport is a perspective projection, which will map the square face of a cube onto a convex quadrilateral. Given two convex quadrilaterals, it is possible to find a planar perspective transformation that maps one to the other. In particular, one can compute a planar perspective map that takes the unit square to an arbitrary convex quadrilateral. This is essentially what you need to map the original image onto the view of the cube face. The image is mapped onto the unit square by dividing image $x$- and $y$-coordinates by its width and height, respectively, as was done above, and then the resulting normalized coordinates mapped to the quadrilateral in the viewport. The coefficients of this planar perspective transformation can be computed by solving a bunch of linear equations, or by pluggin the viewport (projected) coordinates of the cube face into the formulas given here.
