diff --git a/postprocessor.py b/postprocessor.py
index 74a080f..c7c7112 100644
--- a/postprocessor.py
+++ b/postprocessor.py
@@ -865,8 +865,8 @@ class PostProcessor(HDF5Container):
             omega_rho_func, lambda dset: dset["rho"], 1.0 / self._ro.info["unit_time"]
         )
 
-        # Operator to compute the sound speed
-        cs_op = FractionOperator(
+        # Operator to compute the square of the sound speed
+        T_op = FractionOperator(
             lambda dset: dset["P"],
             lambda dset: dset["rho"],
             self._ro.info["unit_velocity"],
@@ -877,21 +877,23 @@ class PostProcessor(HDF5Container):
 
         # Ray tracer for the sound speed
         if self.pp_params.pymses.fft:
-            rt_cs = splatting.SplatterProcessor(
-                self._amr, self._ro.info, cs_op, surf_qty=False
+            rt_T = splatting.SplatterProcessor(
+                self._amr, self._ro.info, T_op, surf_qty=False
             )
         else:
-            rt_cs = raytracing.RayTracer(self._amr, self._ro.info, cs_op)
+            rt_T = raytracing.RayTracer(self._amr, self._ro.info, T_op)
 
         if not self.pp_params.pymses.multiprocessing:
             rt_omega.disable_multiprocessing()
-            rt_cs.disable_multiprocessing()
+            rt_T.disable_multiprocessing()
 
         dmap_omega = rt_omega.process(self._cam[ax_los])
-        dmap_cs = rt_cs.process(self._cam[ax_los])
+        dmap_T = rt_T.process(self._cam[ax_los])
         dmap_col = self.save.root.maps.coldens_z.read()
         map_Q = (
-            (self.lbox * dmap_cs.map.T) * dmap_omega.map.T / (np.pi * self.G * dmap_col)
+            (self.lbox * np.sqrt(dmap_T.map.T))
+            * dmap_omega.map.T
+            / (np.pi * self.G * dmap_col)
         )
         return map_Q