! This PyFerret script calculates the values of section 3.3 of the Paper "Artificial Upwelling - A New Narrative" 
! and was written by M. Juerchott (mjuerchott@geomar.de) on the 12th September 2022.
! In order to run this script without error you have to change the path to the model output files to your personal folder structure.
 
 
! ***  Only run it in -nodisplay -nojnl mode ***
 
 
! Save current state
 
	define REGION/DEFAULT "save"
	set data/save
	set grid/save
 
! Clean everything
 
	cancel data/all
	cancel variable/all
	cancel region/all
 
 
! Script function
 
use /Users/mjuerchott/runs/2021_02_Pipes_noLand_COU85/tavg.02025.01.01.nc               ! [d=1]	ArtUp_8.5
use /Users/mjuerchott/runs/2021_02_noPipes_noLand_COU85_cont/tavg.02025.01.01.nc        ! [d=2]	REF_8.5



! DIC
	let Change_DIC = O_DIC[d=1] - O_DIC[d=2]


! preDIC
	let Change_preDIC = O_preDIC[d=1] - O_preDIC[d=2]


! DICr
	let Change_DICr = O_DICR[d=1] - O_DICR[d=2]

! PP
	let Change_netPP = (O_PHYTNPP[d=1] + O_DIAZNPP[d=1]) - (O_PHYTNPP[d=2] + O_DIAZNPP[d=2])

! ExpP
	let Change_ExpP = O_DETREXP[d=1] - O_DETREXP[d=2]



! Budget calculations

let budget_DIC = Change_DIC[x=@DIN, y=@DIN, z=@DIN, l=16] * 12 * 10^(-15)
list budget_DIC
let budget_DIC_1200 = Change_DIC[x=@DIN, y=@DIN, z=1200:6100@DIN, l=16] * 12 * 10^(-15)
list budget_DIC_1200

let budget_preDIC = Change_preDIC[x=@DIN, y=@DIN, z=@DIN, l=16] * 12 * 10^(-15)
list budget_preDIC
let budget_preDIC_1200 = Change_preDIC[x=@DIN, y=@DIN, z=1200:6100@DIN, l=16] * 12 * 10^(-15)
list budget_preDIC_1200

let budget_remDIC = Change_DICr[x=@DIN, y=@DIN, z=@DIN, l=16] * 12 * 10^(-15)
list budget_remDIC
let budget_remDIC_1200 = Change_DICr[x=@DIN, y=@DIN, z=1200:6100@DIN, l=16] * 12 * 10^(-15)
list budget_remDIC_1200

let budget_CaDIC = budget_DIC - budget_preDIC - budget_remDIC
list budget_CaDIC
let budget_CaDIC_1200 = budget_DIC_1200 - budget_preDIC_1200 - budget_remDIC_1200
list budget_CaDIC_1200



let budget_netPP = Change_netPP[x=@DIN, y=@DIN, z=@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365
list budget_netPP
let budget_ExpP_1200 = Change_ExpP[x=@DIN, y=@DIN, z=1200:6100@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365
list budget_ExpP_1200




let xNPP_REF_85 = O_PHYTNPP[d=2] + O_DIAZNPP[d=2]
let xEXP_REF_85 = O_DETRREMI[d=2]
let xNPP_AU_85 = O_PHYTNPP[d=1] + O_DIAZNPP[d=1]
let xEXP_AU_85 = O_DETRREMI[d=1]

let NPP_REF_85 = xNPP_REF_85[x=@DIN, y=@DIN, z=@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365
let EXP_REF_85 = xEXP_REF_85[x=@DIN, y=@DIN, z=130:6100@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365
let NPP_AU_85 = xNPP_AU_85[x=@DIN, y=@DIN, z=@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365
let EXP_AU_85 = xEXP_AU_85[x=@DIN, y=@DIN, z=130:6100@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365

let xREMI_REF_85 = O_DETRREMI[d=2]
let xREMI_AU_85 = O_DETRREMI[d=1]
let REMI_REF_85_DOWN = xREMI_REF_85[x=@DIN, y=@DIN, z=1200:6100@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365
let REMI_AU_85_DOWN = xREMI_AU_85[x=@DIN, y=@DIN, z=1200:6100@DIN, l=1:16@AVE] * 6.625 * 12 * 10^(-15) *60*60*24*365

! NPP
	list NPP_REF_85
! Export efficiency
	list EXP_REF_85
!Transfer / Sequestration efficiency
	list REMI_REF_85_DOWN
! NPP
	list NPP_AU_85
! Export efficiency
	list EXP_AU_85
!Transfer / Sequestration efficiency
	list REMI_AU_85_DOWN





! Cancel everything
 
	cancel data/all
	cancel variable/all
	cancel region/all
 
! Restore current state
 
	set data/restore
	set grid/restore