When a positron annihilates, two gamma photons are created with orthogonal polarizations. It is possible to use coincidence measurements where both photons undergo Compton scattering to estimate their initial relative polarization orientation. This information is of great interest in gamma imaging systems, such as Positron Emission Tomography, where it may be used as an additional tool to distinguish true coincidence events from scatter and random background. The successful utilization of this principle critically depends on the detector’s angular and energy resolution, which determine its polarimetric performance. In this study, we use Monte Carlo simulations based on the Geant4 toolkit to model two multi-pixel detector configurations identified as prospective for the measurement of gamma-ray polarization in PET. One is based on 2 mm × 2 mm × 20 mm LYSO scintillators and the other is based on 3 mm × 3 mm × 20 mm GAGG scintillators. Each configuration has a pair of modules, each consisting of 64 crystals set up in a single 8 × 8 matrix, where both the recoil electron and the Compton-scattered photon are absorbed. We simulate positron annihilation by generating two back-to-back gamma photons of 511 keV with orthogonal polarizations. The Compton scattering is successfully identified and the modulation of the azimuthal angle difference is clearly observed. The configuration based on GAGG crystals demonstrates slightly better polarimetric performance than the one based on LYSO crystals, reflected in the more pronounced azimuthal modulation.