Ok so I believe the very long way around my answer is below, which is that I think what you’re doing is correct, I just didn’t fully understand it before.
First you pull in your packages and set up a vanilla version of MimiDICE2016:
using Mimi
using MimiDICE2016
m = MimiDICE2016.get_model()
Next you call the Mimi function to set the parameter :TATM to your exogenous temperature
set_param!(m, :TATM, input_temp)
This assigns the values in your input_temp vector to the parameter :TATM in any components with that parameter. In this case the only component is :damages, and you are disconnecting that parameter from the climatedynamics variable :TATM and instead forcing exogenously. Also just note you’ve left the tatm0 parameter alone so it is 0.85 but perhaps that’s fine it’s just your initialization year.
Ok great so at this point you create a marginal model mm.
mm = Mimi.create_marginal_model(m, 1e9)
This basically copies m twice into mm, such that mm contains two models, mm.base and mm.marginal, which at this case point are identical. You’ve also set the delta to be 1e9 so that it knows in the future to divide by that when calculating differences between the two models. By that you are saying your temperature change is as a result of a pulse of 1GtC.
Note here something a bit confusing … MimiDICE has it’s own MimiDICE2016.get_marginal_model which also modifies the emissions. You were smart, don’t use that one, because you don’t want to modify emissions. You are going to exogenously create a different in damages using two different temperature vectors, that is the difference in damages you are looking to see (if I understand correctly).
Anyways, next you said you call the following set_param!(mm, :TATM, new_input_temp) which can’t quite be right since that will error, but here’s what I think you’re doing and what I would do.
update_param!(mm.marginal, :TATM, new_input_temp)
# or set_param! but update_param! is more accurate since there already is one
This will update the :TATM parameter in the :damages component (and any other components with that parameter but there are none here) in the modified sub-model of the mm marginal model.
Now you run the following using the sneaky _compute_scc helper function that the creators didn’t advertise (commented with #“helper function for computing SCC from a MarginalModel, not to be exported or advertised to users”) … but of course you seem to understand what it’s doing so feel free to use it!
scc = MimiDICE2016._compute_scc(mm, year = pulse_year, last_year = 2300, prtp = 0.03, eta = 0.0)
and assuming that the times line up with when your pulse occurred in your exogenous temperature vector is correct, you should be good, except as you said you need to divide by 5 if you want to use this function, because you don’t need to adjust for your pulse being spread out across a timestep. You are correcting for this final line:
scc = sum(df .* marginal_damages * 5) # currently implemented as a 5year step function; so each timestep of discounted marginal damages is multiplied by 5
Or you can write a local copy of the function without the division.
Thus I think you’re doing the right thing, but please tell me if any of this is confusing. Overall I would suggest sanity checking your values, and also perhaps writing your own longer function doing exactly what you want and make sure the values match, just as good practice since I didn’t write up DICE2016 myself nor do I know where your temperature vectors came from 