Merge pull request #124 from brittonsmith/lr

[bugfix] correct ray length for periodic rays

.circleci/config.yml

@@ -11,7 +11,7 @@ commands:
           echo 'TRIDENT_CONFIG=$HOME/.trident/config.tri' >> $BASH_ENV
           echo 'YT_GOLD=8dcc2fa135dfe7a62438507066527caf5cb379c1' >> $BASH_ENV
           echo 'YT_HEAD=yt-4.0' >> $BASH_ENV
-          echo 'TRIDENT_GOLD=test-standard-sph-viz-v5' >> $BASH_ENV
+          echo 'TRIDENT_GOLD=tip' >> $BASH_ENV
           echo 'TRIDENT_HEAD=tip' >> $BASH_ENV
 
   install-dependencies:
@@ -198,14 +198,14 @@ jobs:
 
       - restore_cache:
           name: "Restore answer tests cache."
-          key: test-results-<< parameters.tag >>-<< parameters.yttag >>-v5
+          key: test-results-<< parameters.tag >>-<< parameters.yttag >>-v6
 
       - run-tests:
           generate: 1
 
       - save_cache:
           name: "Save answer tests cache."
-          key: test-results-<< parameters.tag >>-<< parameters.yttag >>-v5
+          key: test-results-<< parameters.tag >>-<< parameters.yttag >>-v6
           paths:
             - ~/answer_test_data/test_results
 

tests/test_light_ray.py

@@ -141,3 +141,13 @@ def test_light_ray_redshift_coverage(self):
         ray = make_simple_ray(ds, start_position=ds.domain_left_edge, end_position=ds.domain_right_edge, lines=['H'])
         assert_almost_equal(ray.r['redshift'][0], 0.00489571, decimal=8)
         assert_almost_equal(ray.r['redshift'][-1], -0.00416831, decimal=8)
+
+    def test_light_ray_redshift_monotonic(self):
+        """
+        Tests to assure a light ray redshift decreases monotonically
+        when ray extends outside the domain.
+        """
+        ds = load(COSMO_PLUS_SINGLE)
+        ray = make_simple_ray(ds, start_position=ds.domain_center,
+                              end_position=ds.domain_center+ds.domain_width)
+        assert((np.diff(ray.data['redshift']) < 0).all())

trident/light_ray.py

@@ -35,7 +35,7 @@
 from yt.data_objects.static_output import Dataset
 
 class LightRay(CosmologySplice):
-    """
+    r"""
     A 1D object representing the path of a light ray passing through a
     simulation.  LightRays can be either simple, where they pass through a
     single dataset, or compound, where they pass through consecutive
@@ -592,6 +592,9 @@ def make_light_ray(self, seed=None, periodic=True,
             for field in all_fields:
                 sub_data[field] = []
 
+            # Keep track of length along full ray.
+            ray_length = ds.quan(0, 'code_length')
+
             # Get data for all subsegments in segment.
             for sub_segment in sub_segments:
                 mylog.info("Getting subsegment: %s to %s." %
@@ -600,12 +603,17 @@ def make_light_ray(self, seed=None, periodic=True,
                 for key, val in field_parameters.items():
                     sub_ray.set_field_parameter(key, val)
                 asort = np.argsort(sub_ray["t"])
-                sub_data['l'].extend(sub_ray['t'][asort] *
-                                     vector_length(sub_ray.start_point,
-                                                   sub_ray.end_point))
-                sub_data['dl'].extend(sub_ray['dts'][asort] *
-                                      vector_length(sub_ray.start_point,
-                                                    sub_ray.end_point))
+                sub_length = vector_length(
+                    sub_ray.start_point,
+                    sub_ray.end_point)
+
+                # redshifts derived from l values
+                sub_data['l'].extend(
+                    sub_ray['t'][asort] * sub_length + ray_length)
+                ray_length += sub_length
+
+                # column densities derived from dl values
+                sub_data['dl'].extend(sub_ray['dts'][asort] * sub_length)
 
                 for field in data_fields:
                     sub_data[field].extend(sub_ray[field][asort])
@@ -660,11 +668,10 @@ def make_light_ray(self, seed=None, periodic=True,
                 sub_data[key] = ds.arr(sub_data[key]).in_cgs()
 
             # Get redshift for each lixel.  Assume linear relation between l
-            # and z.  so z = z_start - (l * (z_range / l_range))
+            # and z.  so z = z_start - z_range * (l / l_range)
             sub_data['redshift'] = my_segment['redshift'] - \
-              (sub_data['l'] * \
-              (my_segment['redshift'] - next_redshift) / \
-              vector_length(my_start, my_end).in_cgs())
+              (sub_data['l'] / ray_length) * \
+              (my_segment['redshift'] - next_redshift)
 
             # When using the peculiar velocity, create effective redshift
             # (redshift_eff) field combining cosmological redshift and