diff --git a/scripts/run.py b/scripts/run.py
index c7b9d0ce06db8d744e17964c75560a0d1a2361ce..489ad8a05a3daeff49059e4f89cd2be46a35ad6d 100644
--- a/scripts/run.py
+++ b/scripts/run.py
@@ -69,7 +69,7 @@ def parse_args():
 	parser.add_argument("--sharpen", default=0, help="Set amount of sharpening applied to NeRF training images. Range 0.0 to 1.0.")
 
 	# parser.add_argument("--tomonerf", action="store_true", help="Apply tomography specific parameters to optimise NeRF output.")
-	parser.add_argument("--tomonerf", type=int, default=0, help="Which transforms to apply [single_channel, exponentiate img data, crop].")
+	parser.add_argument("--tomonerf", type=str, help="Which transforms to apply [single_channel, exponentiate img data, crop].")
 	parser.add_argument("--crop", nargs=3, type=float, help="Percentage to crop unit cube to.")
 
 	return parser.parse_args()
@@ -128,6 +128,13 @@ if __name__ == "__main__":
 
 	if args.tomonerf:
 		print("***Entering NERF for Tomography mode***")
+		print(f"Saving density slices (no crop)")
+		# vertical slices
+		fname = "string"
+		# FIXME - doesn't work from command line yet
+		# testbed.compute_and_save_png_slices()
+		# testbed.compute_and_save_png_slices_vert()
+		# testbed.compute_and_save_png_slices('horiz', flip_y_and_z_axes=True)
 
 		if (format(args.tomonerf, 'b')[-1] == '1'):
 			print("TODO: Single channel mode")
diff --git a/src/python_api.cu b/src/python_api.cu
index d1f4fee279a3ed803c4c5fad140ed3b448eb0c09..299d01d20e5d6137f36f72e46f528afc9d85f4e1 100644
--- a/src/python_api.cu
+++ b/src/python_api.cu
@@ -421,7 +421,7 @@ PYBIND11_MODULE(pyngp, m) {
 		.def("load_file", &Testbed::load_file, py::arg("path"), "Load a file and automatically determine how to handle it. Can be a snapshot, dataset, network config, or camera path.")
 		.def_property("loop_animation", &Testbed::loop_animation, &Testbed::set_loop_animation)
 		.def("compute_and_save_png_slices", &Testbed::compute_and_save_png_slices,
-			py::arg("filename"),
+			py::arg("filename") = "horiz",
 			py::arg("resolution") = Eigen::Vector3i::Constant(256),
 			py::arg("aabb") = BoundingBox{},
 			py::arg("thresh") = std::numeric_limits<float>::max(),
@@ -429,6 +429,15 @@ PYBIND11_MODULE(pyngp, m) {
 			py::arg("flip_y_and_z_axes") = false,
 			"Compute & save a PNG file representing the 3D density or distance field from the current SDF or NeRF model. "
 		)
+		.def("compute_and_save_png_slices_vert", &Testbed::compute_and_save_png_slices,
+			py::arg("filename") = "vert",
+			py::arg("resolution") = Eigen::Vector3i::Constant(256),
+			py::arg("aabb") = BoundingBox{},
+			py::arg("thresh") = std::numeric_limits<float>::max(),
+			py::arg("density_range") = 4.f,
+			py::arg("flip_y_and_z_axes") = true,
+			"Compute & save a PNG file representing the 3D density or distance field from the current SDF or NeRF model. "
+		)
 		.def("compute_and_save_marching_cubes_mesh", &Testbed::compute_and_save_marching_cubes_mesh,
 			py::arg("filename"),
 			py::arg("resolution") = Eigen::Vector3i::Constant(256),
diff --git a/src/testbed_nerf.cu b/src/testbed_nerf.cu
index b0adcfd3e4f2891f78d365d8ec5de9d7ea8939e1..974f28e8ca701dd949b59ff6acfd447a4d145f64 100644
--- a/src/testbed_nerf.cu
+++ b/src/testbed_nerf.cu
@@ -226,14 +226,13 @@ __device__ float network_to_rgb(float val, ENerfActivation activation) {
 
 // No way to modify the derivative for rgb
 __device__ float network_to_rgb_derivative(float val, ENerfActivation activation) {
-	return 0.0f;
-	// switch (activation) {
-	// 	case ENerfActivation::None: return 1.0f;
-	// 	case ENerfActivation::ReLU: return val > 0.0f ? 1.0f : 0.0f;
-	// 	case ENerfActivation::Logistic: { float density = tcnn::logistic(val); return density * (1 - density); };
-	// 	case ENerfActivation::Exponential: return __expf(tcnn::clamp(val, -10.0f, 10.0f));
-	// 	default: assert(false);
-	// }
+	switch (activation) {
+		case ENerfActivation::None: return 1.0f;
+		case ENerfActivation::ReLU: return val > 0.0f ? 1.0f : 0.0f;
+		case ENerfActivation::Logistic: { float density = tcnn::logistic(val); return density * (1 - density); };
+		case ENerfActivation::Exponential: return __expf(tcnn::clamp(val, -10.0f, 10.0f));
+		default: assert(false);
+	}
 }
 
 __device__ float network_to_density(float val, ENerfActivation activation) {
@@ -261,9 +260,9 @@ __device__ float network_to_density_derivative(float val, ENerfActivation activa
 // Ignore neurons 0, 1 and 2!
 __device__ Array3f network_to_rgb(const tcnn::vector_t<tcnn::network_precision_t, 4>& local_network_output, ENerfActivation activation) {
 	return {
-		network_to_rgb(0.0f, activation),
-		network_to_rgb(0.0f, activation),
-		network_to_rgb(0.0f, activation)
+        network_to_rgb(float(local_network_output[0]), activation),
+        network_to_rgb(float(local_network_output[1]), activation),
+        network_to_rgb(float(local_network_output[2]), activation)
 	};
 }