Add Generation of Restart script

util/gener_perturb_restart.bash

@@ -0,0 +1,85 @@
+#!/bin/bash
+#
+# -- Author : François Massonnet, francois.massonnet@ic3.cat
+# -- Date   : 30 Jan 2015
+# -- At     : IC3, Barcelona
+# -- Modified : 19 Jan 2016, omar.bellprat@bsc.es 
+#
+# -- Purpose: Generation of an arbitrary number of NEMO oceanic restarts that are copies of a reference, plus a perturbation
+#
+# -- Method : The reference file is duplicated in this script, then read by a R script. The perturbation is introduced and finally written to this latter file.
+#             The script must be placed into the restart directory (e.g. NEMO_Restart_23). The generated restarts have to be renamed after generation. This 
+#             script has been tested on MareNostrum3.
+#
+# -- Input  : NEMO ocean restart from an EC-Earth 3.1 run
+# -- Output : N restarts with the same name,but with an index fc0, fc1, ... fcN-1 appended
+#
+# -- Limitations: Only the surface conditions are perturbed (level index: 1) but this can be changed in the R script
+
+module load R
+
+set -o errexit
+set -o nounset
+set -x
+
+
+if [ $# == 0 ] ; then
+  echo "gener_perturbation.bash ocean_restart_file Nmembers"
+  exit
+fi
+
+filein=$1
+nmemb=$2
+
+# ---------------------------------------------------------
+
+var=tn          # Variable to be perturbed
+per=0.0001      # Standard deviation of gaussian perturbation to be applied,
+                # in units of the variable (for tn: in K for example)
+
+cp $filein ${filein}.backup # Do a backup
+
+for jmemb in `seq 0 $(( $nmemb -1 ))`
+do
+  echo $jmemb
+  # 1. Make a copy of the original file, with the new name
+  filenew="${filein%.nc}_fc${jmemb}.nc"
+  cp $filein $filenew
+
+  # 2. Prepare the R script
+
+  echo "#!/usr/bin/env Rscript
+  library(ncdf)
+
+  # François Massonnet, 30 Jan 2015
+  # Adds a gaussian perturbation at the first level of a 3D field
+  # Tested only for NEMO restarts
+  #
+  # This script should be called by a bash script so that the variable and file names are specified, as well as the perturbation
+
+  varname='$var'
+  filein <- '$filenew'
+  ex.nc           <- open.ncdf(filein,write=TRUE)
+  spert <- $per
+
+  myvar     <- get.var.ncdf(ex.nc, varname)
+  myvarpert <- myvar
+  for (i in seq(1,dim(myvar)[1])){
+    for (j in seq(1,dim(myvar)[2])){
+      if (myvar[i,j,1] != 0){
+        myvarpert[i,j,1] = myvarpert[i,j,1] + rnorm(1,sd=spert)
+      }
+    }
+  }
+
+  put.var.ncdf(ex.nc,varname,myvarpert)
+  close.ncdf(ex.nc)" > Rtmpfile.R
+
+  chmod 744 Rtmpfile.R
+
+  # 3. Run the R script, that produces the new NetCDF
+  ./Rtmpfile.R 
+
+  rm -f Rtmpfile.R
+done
+