Recently I've been trying to implement OpenGL perspective projection matrix, but stuck at understanding of how to derive projected Z value. Here's an example of it
Derivation of $x$ and $y$ values is pretty straightforward. Firstly, I need to project point onto near plane like so: $x_p=nx/z$ where $n$ is distance to near plane, and remap $x$ from $l<=x<=r$ to $-1<=x<=1$. Same steps for $y$, but $z$ has differences. I cannot just perform same steps for it as I did for $x$ and $y$ by solving simple inequalities. Question here: why I cannot just simply calculate it by solving inequality of type $n<=z<=f$, but rather must solve these equations: $-n=pn+q$ and $f=pf+q$ (assuming I want to map $z$ into $[-1;1]$ range)? Is there some geometric meaning for that?
Source for matrix derivation that I was using: https://www.codeguru.com/cpp/misc/misc/graphics/article.php/c10123/Deriving-Projection-Matrices.htm#page-3

I guess I've found the answer for this question. The problem is that I missed a point of normalization and clipping. Before implementing frustum view volume normalization it's a good idea to map frustum view volume to orthographic view volume and having that done map one to normalized cube which is NDC. That's why in all books firstly we derive orthographic projection matrix and then perspective projection.
And about "condensing" $z$ between near and far planes. It's pretty straightforward. In order to obtain values for orthographic view volume we need 1) normalize view frustum range 2) stretch normalized range to orthographic view volume of form $[0, far]$. To do so first we need to offset our segment $[near, far]$ to the origin. Doing that we get $[0, far-near]$. The only thing that is left is to normalize the latest component which is simple division by $far-near$. Thus we get normalized coordinates. Now we can easily obtain orthographic view volume range by just multiply everything what we have gotten by $far$ value which is $[0, far]$ and it's how we get $\frac {f}{f-n} - \frac {fn}{f-n}$