openalea / WALTer

WALTer is a 3D FSPM Wheat model
Other
4 stars 8 forks source link

Scene and light units #38

Closed christian34 closed 6 years ago

christian34 commented 6 years ago

Switch walter to consistent input unit for caribu: output / light responses should be checked before merging

Thris branch was based on PR #17: it may be good to merge that branch first

Emblanc commented 6 years ago

When running WALTer with this PR, the values for the intercepted PAR are lower than before (by a 10⁴ factor), which is consistent with a change of units. On the face of it, these changes seem OK to me since we knew our plants used to intercept too much light (by approximatively a 10⁴ factor). But I will need to confirm that these changes of values are indeed correct biologically speaking.

This PR passes all the usual tests on my computer. However, when I look at the intercepted PAR, there are very important changes between one day and the next day for some plants (seems like the issue #14 is back). My guess is that the bug is present but that it was not detected by our test (see #40)

christian34 commented 6 years ago

I extend the test to all axes, normalising for incident par : I found mean variation less than 50% for all axes: @Emblanc : is it okay ? We can also refine the test to check min / max variation range

Emblanc commented 6 years ago

The new test seems great ! Unfortunately, the PR still passes all the usual tests on my computer, which is not consistent with what I am observing when running the test simulation and looking at the PAR_per_axes.csv output file myself...

Emblanc commented 6 years ago

I took a closer look at it and it seems that the bug has not been detected because it happens on day 204 - 205 and the test is only conducted on a few days : it is not conducted on these specific days.

christian34 commented 6 years ago

Okay cool that you locate the bug, I will have a look at these specific days !

Emblanc commented 6 years ago

Here is a document showing the results I have when I run the sim_scheme_test. The results are the same as in the test repository. Shift_in_light.pdf

christian34 commented 6 years ago

Okay, I investigate a bit more on what happen on day 205 for tiller (1,2,1) of plant 3. The high PAR values comes from the way sumPAR is computed: 1) on relevant organ, PAR is computed as : Ei (PAR irradiance) * area_of_the_primitive 2) tiller_surface is computed as the sum of visible_area of all relevant organ belonging to an axis (with math formulae + estimate of visibility) 3) organ PAR are summed per axis and divided by tiller_surface and Temperature ???, to get a 'corrected for visibility (and temperature ??)' irradiance (called sum_PAR)

The strange values for sum_PAR occur when the tiller_surface << sum(primitive_area). For example, at the considered time step, the tiller_surface of axe (1,2,1) of plant 3 is less than 1 % of the sum of primitive organ area, whereas for other axes it is from 3 to 50 %.

These discrepancies can come from discretisation effects (math cylinder area are not equal to primitive area) and bugs/ problems in visibility computation (it is difficult to assess a priori what part of the primitive is exposed to light when geometric object are intricated and when light comes from several directions)

Probably the best solution will be to switch to the representation of 'only visible' part of organs (PR #9) and compute tiller area as the sum of contributing primitives

@Emblanc Meanwhile, we can try to limit the bug by limiting the correction to cases where tiller_surface is at least a few % of primitive area (let say 5%) (and then work on PR #9)?

Whatever the solution kept, I will have to convert tiller_surface to square meter before dividing organ PAR

christian34 commented 6 years ago

Just note that I updated my previous comment due to a mistake in result aggregation : the problem with tiller_surface is by far less dramatic than what my original comment suggested...

christian34 commented 6 years ago

Using correct unit for tiller_surface (cm2 -> m2) in tiller irradiance computation allows to pass the test of maximal variation < 100 for all plants for the whole simulation (the problematic tiller at step 205 disappear, due to proper selection of surviving axis after unit correction). I suspect that shift in light bug can re-appear for other simulation, but this can be discuss longer in PR #9.

To me the PR can be merged: @Emblanc @pradal is it okay for you ?

Emblanc commented 6 years ago

I am checking light outputs to see if everything is OK : I saved the intercepted PAR for each axis on each day and I compared it with the incident PAR. It seems that the intercepted PAR is sometimes higher than the incident PAR, which is abnormal. graph_intercepted_par But it is likely that I am doing something wrong when I check the outputs. Here is how I did it : 1) I saved the intercepted PAR for each axis and each day (same as line 2348 in WALTer.lpy but not divided by temperature or by surface ; so unit = µmolPAR) 2) I compared it with currentPAR(µmol/m²)*PlotArea(m²) -> unit = µmolPAR

christian34 commented 6 years ago

I add some lines to automate this check (see check_light_balance in test_light_interception), and found, like you, that total interception is higher than incident PAR on plot_area. This situation may however be perfectly normal : 1) by convention, incoming light flux is quantified as the light flux density passing through a horizontal surface (ie only the vertical flux is accounted). The total light flux emitted by sources is higher than the vertical component, due to cosine effects (for example an almost horizontal light source will have very low vertical component, whatever its intensity). It is therefore possible that total interception > total vertical flux 2) the reference simulation is run with infinite = False, which increase effect of 1) on border plants (more 'horizontal light' is passing through an isolated canopy) and allows ray lights falling 'outside' the soil below the plot to be intercepted by plants

Generally, we can observe (but some extreme case can exists) that total interception < incoming vertical light on the plot when infinite is True (and this is effectively the case if you run check_light_balance with infinity_CARIBU=1). In this situation, you will also always guarantee that vertical flux through the soil is below incoming light (but to test that, you will have to add primitive representing the soil).

Emblanc commented 6 years ago

Thank you for this complete explanation. I'm glad that everything seems normal after all. I think this PR is ready to be merged then.