We discuss peculiarities that arise in the computation of real-emission contributions to observables that contain Heaviside functions. A prominent example of such a case is the zero-jettiness soft function in SCET, whose calculation at next-to-next-to-next-to-leading order in perturbative QCD is an interesting problem. Since the zero-jettiness soft function distinguishes between emissions into different hemispheres, its definition involves θ-functions of light-cone components of emitted soft partons. This prevents a direct use of multi-loop methods, based on reverse unitarity, for computing the zero-jettiness soft function in high orders of perturbation theory. We propose a way to bypass this problem and illustrate its effectiveness by computing various non-trivial contributions to the zero-jettiness soft function at NNLO and N3LO in perturbative QCD.