The Snow Grain Size and Pollution (SGSP) algorithm is a method applied to Moderate Resolution Imaging Spectroradiometer data to estimate snow grain size from space-borne measurements. Previous studies validate and quantify potential sources of error in this method, but because it assumes flat snow surfaces, however, large scale variations in surface normals can cause biases in its estimates due to its dependence on solar and observation zenith angles. To address these variations, we apply the Monte Carlo method for photon transport using data containing the single scattering properties of different ice crystals to calculate polarization states of reflected monochromatic light at 1500nm from modeled snow surfaces. We evaluate the dependence of these polarization states on solar and observation geometry at 1500nm because multiple scattering is generally a mechanism for depolarization and the ice crystals are relatively absorptive at this wavelength. Using 1500nm thus results in a higher number of reflected photons undergoing fewer scattering events, increasing the likelihood of reflected light having higher degrees of polarization. In evaluating the validity of the model, we find agreement with previous studies pertaining to near-infrared spectral directional hemispherical reflectance (i.e. black-sky albedo) and similarities in measured bidirectional reflectance factors, but few studies exist modeling polarization states of reflected light from snow surfaces. Here, we present novel results pertaining to calculated polarization states and compare dependences on solar and observation geometry for different idealized snow surfaces. If these dependencies are consistent across different ice particle shapes and sizes, then these findings could inform the SGSP algorithm by providing useful relationships between measurable physical quantities and solar and observation geometry to better understand variations in snow surface normals from remote sensing observations.