source: TOOLS/CMIP6_FORCING/SOLAR/process.sh @ 3381

Last change on this file since 3381 was 3381, checked in by jgipsl, 4 years ago

Initial version of scripts to create solar forcings.

O Boucher

  • Property svn:executable set to *
File size: 7.0 KB
Line 
1#--this is a script that does it all
2
3dirout="OUTPUT"
4
5#--first clean up things
6rm -f solarforcing_picontrol_fx_3.1_SSI_regridded_2_TSI_anomaly.nc
7rm -f solarforcing_picontrol_fx_3.1_SSI_regridded_6_TSI_anomaly.nc
8rm -f solarforcing_ref_day_3.1_SSI_regridded_2_TSI_anomaly.nc
9rm -f solarforcing_ref_day_3.1_SSI_regridded_6_TSI_anomaly.nc
10rm -f solarforcing_ext_day_3.1_SSI_regridded_2_TSI_anomaly.nc
11rm -f solarforcing_ext_day_3.1_SSI_regridded_6_TSI_anomaly.nc
12
13#--run matlab script for piControl
14matlab -r "run  CMIP6_piControl_2bands.m, quit"
15matlab -r "run  CMIP6_piControl_6bands.m, quit"
16
17#--run matlab script for Historical and Future
18matlab -r "run  CMIP6_historical_2bands.m, quit"
19matlab -r "run  CMIP6_historical_6bands.m, quit"
20
21#--run matlab script for Future Sensitivity
22matlab -r "run  CMIP6_sensitivity_2bands.m, quit"
23matlab -r "run  CMIP6_sensitivity_6bands.m, quit"
24
25output=${dirout}"/Control"
26if [ ! -d "$output" ]; then
27mkdir $output
28fi
29
30#--as computed by CMIP6 for piControl
31filein2bands=./solarforcing_picontrol_fx_3.1_SSI_regridded_2_TSI_anomaly.nc
32filein6bands=./solarforcing_picontrol_fx_3.1_SSI_regridded_6_TSI_anomaly.nc
33
34cdo duplicate,360 $filein2bands ${output}/solarforcing_v3.1_piControl_ave_2bands_360days.nc
35cdo duplicate,360 $filein6bands ${output}/solarforcing_v3.1_piControl_ave_6bands_360days.nc
36cdo duplicate,365 $filein2bands ${output}/solarforcing_v3.1_piControl_ave_2bands_365days.nc
37cdo duplicate,365 $filein6bands ${output}/solarforcing_v3.1_piControl_ave_6bands_365days.nc
38cdo duplicate,366 $filein2bands ${output}/solarforcing_v3.1_piControl_ave_2bands_366days.nc
39cdo duplicate,366 $filein6bands ${output}/solarforcing_v3.1_piControl_ave_6bands_366days.nc
40
41#--now the same duplication but to match Gregorian calendar
42outputG=${dirout}"/ControlGregorian"
43if [ ! -d "$outputG" ]; then
44mkdir $outputG
45fi
46
47for year in {1850..2500}
48do
49nbdays=`python -c "import calendar ; print 366 if calendar.isleap(${year}) else 365"`
50rm -f ${outputG}/solarforcing_v3.1_piControl_ave_2bands_${year}.nc
51rm -f ${outputG}/solarforcing_v3.1_piControl_ave_6bands_${year}.nc
52ln -s ../Control/solarforcing_v3.1_piControl_ave_2bands_${nbdays}days.nc ${outputG}/solarforcing_v3.1_piControl_ave_2bands_${year}.nc
53ln -s ../Control/solarforcing_v3.1_piControl_ave_6bands_${nbdays}days.nc ${outputG}/solarforcing_v3.1_piControl_ave_6bands_${year}.nc
54done
55
56rm -f ${filein2bands} ${filein6bands}
57
58#--present-day average for unofficial CMIP6 pdControl_ave
59
60filein2bands=./solarforcing_ref_day_3.1_SSI_regridded_2_TSI_anomaly.nc
61filein6bands=./solarforcing_ref_day_3.1_SSI_regridded_6_TSI_anomaly.nc
62
63#year1=1995
64#year2=2004
65year1=1991
66year2=2010
67start_time=`python -c "from datetime import date; print (date($year1,1,1)-date(1850,1,1)).days+1"`
68end_time=`python -c "from datetime import date; print (date($year2,12,31)-date(1850,1,1)).days"`
69echo $year1 $year2 $start_time $end_time
70
71fileout2bands=solarforcing_ref_day_3.1_SSI_regridded_2_TSI_anomaly_${year1}-${year2}.nc
72fileout6bands=solarforcing_ref_day_3.1_SSI_regridded_6_TSI_anomaly_${year1}-${year2}.nc
73
74#--somewhat dirty average as in principle ssi_frac should be weighted by tsi !
75#--note that seltimestep starts with 1
76cdo timmean -seltimestep,$start_time,$end_time $filein2bands $fileout2bands
77cdo timmean -seltimestep,$start_time,$end_time $filein6bands $fileout6bands
78
79#--duplicate 360 365 or 366 values for annual file
80cdo duplicate,360 $fileout2bands ${output}/solarforcing_v3.1_pdControl_ave_2bands_360days.nc
81cdo duplicate,360 $fileout6bands ${output}/solarforcing_v3.1_pdControl_ave_6bands_360days.nc
82cdo duplicate,365 $fileout2bands ${output}/solarforcing_v3.1_pdControl_ave_2bands_365days.nc
83cdo duplicate,365 $fileout6bands ${output}/solarforcing_v3.1_pdControl_ave_6bands_365days.nc
84cdo duplicate,366 $fileout2bands ${output}/solarforcing_v3.1_pdControl_ave_2bands_366days.nc
85cdo duplicate,366 $fileout6bands ${output}/solarforcing_v3.1_pdControl_ave_6bands_366days.nc
86
87#--now the same duplication but to match Gregorian calendar
88for year in {1850..2500}
89do
90nbdays=`python -c "import calendar ; print 366 if calendar.isleap(${year}) else 365"`
91rm -f ${outputG}/solarforcing_v3.1_pdControl_ave_2bands_${year}.nc
92rm -f ${outputG}/solarforcing_v3.1_pdControl_ave_6bands_${year}.nc
93ln -s ../Control/solarforcing_v3.1_pdControl_ave_2bands_${nbdays}days.nc ${outputG}/solarforcing_v3.1_pdControl_ave_2bands_${year}.nc
94ln -s ../Control/solarforcing_v3.1_pdControl_ave_6bands_${nbdays}days.nc ${outputG}/solarforcing_v3.1_pdControl_ave_6bands_${year}.nc
95done
96
97rm -f $fileout2bands $fileout6bands
98
99#--split files for historical
100filein2bands=./solarforcing_ref_day_3.1_SSI_regridded_2_TSI_anomaly.nc
101filein6bands=./solarforcing_ref_day_3.1_SSI_regridded_6_TSI_anomaly.nc
102
103output=${dirout}"/Historical"
104if [ ! -d "$output" ]; then
105mkdir $output
106fi
107
108#--note that the time axis starts from 0 with ncks.... so I put -1 here
109end_time=-1
110
111for year in {1850..2014}
112do
113yearp1=$((year+1))
114nbday=`python -c "from datetime import date; print (date($yearp1,1,1)-date($year,1,1)).days"`
115start_time=$((end_time+1))
116end_time=$((start_time+nbday-1))
117echo $year $yearp1 $start_time $end_time
118rm -f ${output}/solarforcing_v3.1_daily_2bands_${year}.nc
119rm -f ${output}/solarforcing_v3.1_daily_6bands_${year}.nc
120ncks -d time,$start_time,$end_time $filein2bands ${output}/solarforcing_v3.1_daily_2bands_${year}.nc
121ncks -d time,$start_time,$end_time $filein6bands ${output}/solarforcing_v3.1_daily_6bands_${year}.nc
122done 
123
124#--split files for Future
125output=${dirout}"/Future"
126if [ ! -d "$output" ]; then
127mkdir $output
128fi
129
130for year in {2015..2299}
131do
132yearp1=$((year+1))
133nbday=`python -c "from datetime import date; print (date($yearp1,1,1)-date($year,1,1)).days"`
134start_time=$((end_time+1))
135end_time=$((start_time+nbday-1))
136echo $year $yearp1 $start_time $end_time
137rm -f ${output}/solarforcing_v3.1_daily_2bands_${year}.nc
138rm -f ${output}/solarforcing_v3.1_daily_6bands_${year}.nc
139ncks -d time,$start_time,$end_time $filein2bands ${output}/solarforcing_v3.1_daily_2bands_${year}.nc
140ncks -d time,$start_time,$end_time $filein6bands ${output}/solarforcing_v3.1_daily_6bands_${year}.nc
141done
142
143rm -f $filein2bands $filein6bands
144
145#--split files for Future Sensitivity
146filein2bands=./solarforcing_ext_day_3.1_SSI_regridded_2_TSI_anomaly.nc
147filein6bands=./solarforcing_ext_day_3.1_SSI_regridded_6_TSI_anomaly.nc
148
149output=${dirout}"/FutureSensitivity"
150if [ ! -d "$output" ]; then
151mkdir $output
152fi
153
154#--reinitialising end_time to 31/12/2014
155end_time=`python -c "from datetime import date; print (date(2014,12,31)-date(1850,1,1)).days"`
156
157for year in {2015..2299}
158do
159yearp1=$((year+1))
160nbday=`python -c "from datetime import date; print (date($yearp1,1,1)-date($year,1,1)).days"`
161start_time=$((end_time+1))
162end_time=$((start_time+nbday-1))
163echo $year $yearp1 $start_time $end_time
164rm -f ${output}/solarforcing_v3.1_daily_2bands_${year}.nc
165rm -f ${output}/solarforcing_v3.1_daily_6bands_${year}.nc
166ncks -d time,$start_time,$end_time $filein2bands ${output}/solarforcing_v3.1_daily_2bands_${year}.nc
167ncks -d time,$start_time,$end_time $filein6bands ${output}/solarforcing_v3.1_daily_6bands_${year}.nc
168done
169
170rm -f $filein2bands $filein6bands
Note: See TracBrowser for help on using the repository browser.