In the context of Multiscale Hybrid-Mixed (MHM) methods for linear elasticity, an attractive strategy for solving the local Neumann subproblems is to employ a Galerkin method on continuous piecewise polynomial spaces for their primal formulation. While computationally efficient, this approach generally yields approximate stress fields that are not $H(div,\Omega)$-conforming, which compromises local conservation and limits the quality of the stress approximation. In the simulation of coupled geomechanics and fluid flow problems in porous media, such a limitation may significantly affect the overall quality and reliability of the numerical results. To address this limitation, we develop $H(div,\Omega)$-conforming post-processing techniques for MHM methods with both weakly and strongly enforced symmetry of the stress variable. The MHM framework decomposes the exact solution into global and local contributions, leading after discretization to a global skeletal formulation coupled with independent local Neumann subproblems. We introduce two element-wise post-processing reconstruction strategies defined on submeshes, designed to recover $H(div,\Omega)$-conforming stress tensors while enforcing either weak or strong symmetry. We prove that the reconstructed stress tensors converge optimally in the $H(div,\Omega)$-norm with respect to the diameters of the face mesh and the local submeshes. Moreover, we show that the local divergence projection coincides with the piecewise continuous polynomial projection of the source term onto the submesh, yielding optimal convergence in the $L^2(\Omega)$-norm with respect to the local mesh size. Numerical experiments confirm the theoretical results and demonstrate the robustness and effectiveness of the proposed strategies in multilayer elasticity problems, aimed at real-world applications such as faulted subsurface reservoirs.