A more general strategy based on the application of Wigner-Eckart's theorem to transition density has been outlined by McWeeny . Here we follow McWeeny's proposal and develop second-quantization formulation of this approach. The essential components of our formalism are very similar to the implementation of SOC calculations within the restricted active space state interaction (RAS-SI) framework in the Molcas program . In this approach, one needs to compute only one multiplet component (the M s =0 state) from which the reduced spinless density matrix is generated. This density matrix is used to compute a reduced matrix element and the entire SOC matrix is obtained by the application of Wigner-Eckart's theorem. The phase problem is avoided by construction, because the density matrix corresponding to a reduced matrix element is computed from only one transition between the pair of target states (or between the reference and the target EOM state). The theory is illustrated by application to the EOM-CCSD wave functions. As numerical exam-ples, we report calculations for a set of small open-shell molecules, considered as candidates for laser cooling experiments, a set of diradicals used as benchmark systems in previous studies , several transition-metal ions, and a Fe(II) SMM with a large spin-reversal barrier (magnetic anisotropy) for which a synthetic analogue has been synthesized .