! This PyFerret script plots FIG_S1 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_11_Pipes_noLand_COU60/tavg.02025.01.01.nc 	![d=2]	ArtUp_6.0
	use /Users/mjuerchott/runs/2021_11_Pipes_noLand_COU45/tavg.02025.01.01.nc 	![d=3]	ArtUp_4.5
	use /Users/mjuerchott/runs/2021_11_Pipes_noLand_COU26/tavg.02025.01.01.nc  	![d=4]	ArtUp_2.6
	use /Users/mjuerchott/runs/2021_02_Pipes_noLand_CTRL/tavg.02025.01.01.nc 	![d=5]	ArtUp_0.0


! pipe area over time

	let pipe_mask_85 = if O_ZPIPE[d=1] GT 50 then 1 else 0
	let pipe_mask_60 = if O_ZPIPE[d=2] GT 50 then 1 else 0
	let pipe_mask_45 = if O_ZPIPE[d=3] GT 50 then 1 else 0
	let pipe_mask_26 = if O_ZPIPE[d=4] GT 50 then 1 else 0
	let pipe_mask_00 = if O_ZPIPE[d=5] GT 50 then 1 else 0

	let pipe_area_m2_85 = (G_AREAT[d=1] * pipe_mask_85)
	let pipe_area_m2_60 = (G_AREAT[d=2] * pipe_mask_60)
	let pipe_area_m2_45 = (G_AREAT[d=3] * pipe_mask_45)
	let pipe_area_m2_26 = (G_AREAT[d=4] * pipe_mask_26)
	let pipe_area_m2_00 = (G_AREAT[d=5] * pipe_mask_00)

	let/unit="[km<&sup2>]" pipe_area_km2_85 = pipe_area_m2_85[i=@SUM,j=@SUM] / 1000000
	let/unit="[km<&sup2>]" pipe_area_km2_60 = pipe_area_m2_60[i=@SUM,j=@SUM] / 1000000
	let/unit="[km<&sup2>]" pipe_area_km2_45 = pipe_area_m2_45[i=@SUM,j=@SUM] / 1000000
	let/unit="[km<&sup2>]" pipe_area_km2_26 = pipe_area_m2_26[i=@SUM,j=@SUM] / 1000000
	let/unit="[km<&sup2>]" pipe_area_km2_00 = pipe_area_m2_00[i=@SUM,j=@SUM] / 1000000



! pipe length average over time

	let pipe_length_2D_85 = O_ZPIPE[d=1] * pipe_mask_85
	let pipe_length_2D_60 = O_ZPIPE[d=2] * pipe_mask_60
	let pipe_length_2D_45 = O_ZPIPE[d=3] * pipe_mask_45
	let pipe_length_2D_26 = O_ZPIPE[d=4] * pipe_mask_26
	let pipe_length_2D_00 = O_ZPIPE[d=5] * pipe_mask_00

	let/unit="[m]" pipe_lenght_ave_85 = pipe_length_2D_85[i=@DIN,j=@DIN] / pipe_mask_85[i=@DIN,j=@DIN]
	let/unit="[m]" pipe_lenght_ave_60 = pipe_length_2D_60[i=@DIN,j=@DIN] / pipe_mask_60[i=@DIN,j=@DIN]
	let/unit="[m]" pipe_lenght_ave_45 = pipe_length_2D_45[i=@DIN,j=@DIN] / pipe_mask_45[i=@DIN,j=@DIN]
	let/unit="[m]" pipe_lenght_ave_26 = pipe_length_2D_26[i=@DIN,j=@DIN] / pipe_mask_26[i=@DIN,j=@DIN]
	let/unit="[m]" pipe_lenght_ave_00 = pipe_length_2D_00[i=@DIN,j=@DIN] / pipe_mask_00[i=@DIN,j=@DIN]



! plot routine

	cancel windows/all
	set windows/new
	set text/font='Times New Roman'
	set text/isiz=4

	go landscape2x2.jnl
	define view/xlimits=0.2,0.7/ylimits=0.05,0.55 ll4
	define view/xlimits=0.2,0.7/ylimits=0.5,1 ul4


	set view ul4
	plot/nolabel/VLIMIT=100000000:300000000:50000000/HLIMITS=2020:2101:10/color=(100,68,0)/thick=2 pipe_area_km2_00
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=green/thick=2 pipe_area_km2_26
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=red/thick=2 pipe_area_km2_45
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=blue/thick=2 pipe_area_km2_60
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=black/thick=2 pipe_area_km2_85
	label/nouser -0.8,-0.8,0,0,0.2 "(a)"
	label/nouser -1.3,1.7,0,90,0.13 "Area [km<&sup2>]"
	label/nouser 2.5,-0.8,0,0,0.13 "[years]"


	set view ll4
	plot/nolabel/VLIMIT=400:1100:100/HLIMITS=2020:2101:10/color=(100,68,0)/thick=2 pipe_lenght_ave_00
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=green/thick=2 pipe_lenght_ave_26
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=red/thick=2 pipe_lenght_ave_45
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=blue/thick=2 pipe_lenght_ave_60
	plot/nolabel/overlay/HLIMITS=2020:2101:10/color=black/thick=2 pipe_lenght_ave_85
	label/nouser -0.8,-0.8,0,0,0.2 "(b)"
	label/nouser -1.3,1.7,0,90,0.13 "Pipe Source Depth [m]"
	label/nouser 2.5,-0.8,0,0,0.13 "[years]"


! plot legend

	set text/font='Times New Roman'/color=black
	label/nouser 0.2,7.9,-1,0,0.11 - RCP 8.5

	set text/font='Times New Roman'/color=blue
	label/nouser 0.2,7.6,-1,0,0.11 - RCP 6.0

	set text/font='Times New Roman'/color=red
	label/nouser 0.2,7.3,-1,0,0.11 - RCP 4.5

	set text/font='Times New Roman'/color=green
	label/nouser 0.2,7.0,-1,0,0.11 - RCP 2.6

	set text/font='Times New Roman'/color=(100,68,0)
	label/nouser 0.2,6.7,-1,0,0.11 - no emission

	set text/font='Times New Roman'/color=black



frame/file = "FIG_S1.png"



! Cancel everything

        cancel data/all
        cancel variable/all
        cancel region/all

! Restore current state

        set data/restore
        set grid/restore