Geomaterials exhibit complex failure mechanisms characterized by strain localization and discontinuities (such as micro-crack formation and propagation), posing substantial challenges in their numerical modelling using the continuum-based methods. These limitations are typically addressed through algorithmic interventions or by using a non-local formulation. Peridynamic is one such method that inherently overcomes these limitations by replacing the partial differentials with non-local integral equations, enabling material points to interact with neighbouring points within a defined horizon. Similar to other numerical approaches, the application of peridynamic in geomechanics necessitates precise calibration of elastic parameters, as they play a crucial role in governing the plastic behavior of geomaterials. In this context, the present investigation studied the optimization of non-ordinary state-based peridynamic formulations using a geomaterial test specimen with a 1:2 aspect ratio under compression to evaluate the influence of critical numerical parameters, namely horizon size and material points discretization, on the accuracy of predicted elastic modulus. The results highlighted the necessity of selecting optimal combinations of mesh density and horizon size to achieve convergence toward input elastic properties. Furthermore, the obtained optimized parameters were used to simulate a series of plane strain compression tests on geomaterials to gain insight into plastic deformation and shear band formation. The study affirmed that parameter calibration is fundamental for accurately capturing both elastic and plastic behaviors of geomaterials. This calibrated model can offer significant potential for modelling failure surfaces below foundations, behind a retaining wall, and on a slope.