Difference between revisions of "Tool Box"

./makelmdz_fcm -arch my_arch_file -p std -d 64x48x30 newstart

./makelmdz_fcm -arch my_arch_file -p std -d 64x48x30 start2archive

from    numpy import *
import  numpy                 as        np
import  matplotlib.pyplot     as        mpl
import  math
import xarray as xr

#FIRST WE GET THE DATA FROM THE GCM SIMULATION
nc = xr.open_dataset('startfi.nc',decode_times=False) # can be any netcdf file (e.g. start/startfi.nc files)

# BELOW PHYSICAL VARIABLES
physical_points=nc['physical_points']
lat=nc['latitude']
lon=nc['longitude']
aire_GCM=nc['area']

# BELOW THE VARIABLE WE WANT TO UPDATE
tsurf=nc['tsurf']
new_tsurf = np.empty(len(physical_points))

# LOOP TO MODIFY THE VARIABLE
for i in range(0,len(physical_points),1):
    new_tsurf[i]=300. # here you put whatever you want ; in this exemple, we assume an isothermal temperature distribution

nc['tsurf'].values = new_tsurf
nc.to_netcdf('restartfi.nc')

# SANITY CHECK PLOTS

fig = mpl.figure(1)
mpl.plot(physical_points,tsurf)
mpl.plot(physical_points,new_tsurf)
mpl.xlabel('GCM Physical Points')
mpl.ylabel('Tsurf (K)')
mpl.show()

trunk/LMDZ.GENERIC/utilities

./streamfunction.e

mv restart.nc start.nc
mv restartfi.nc startfi.nc

#!/bin/bash
###########################################################################
# Script to perform several chained LMD Mars GCM simulations
# SET HERE the maximum total number of simulations

nummax=100

###########################################################################


echo "---------------------------------------------------------"
echo "STARTING LOOP RUN"
echo "---------------------------------------------------------"

dir=`pwd`
machine=`hostname`
address=`whoami`

# Look for file "num_run" which should contain 
# the value of the previously computed season
# (defaults to 0 if file "num_run" does not exist)
if [[ -r num_run ]] ; then
  echo "found file num_run"
  numold=`cat num_run`
else
  numold=0
fi
echo "numold is set to" ${numold}


# Set value of current season 
(( numnew = ${numold} + 1 ))
echo "numnew is set to" ${numnew}

# Look for initialization data files (exit if none found)
if [[ ( -r start${numold}.nc  &&  -r startfi${numold}.nc ) ]] ; then
   \cp -f start${numold}.nc start.nc
   \cp -f startfi${numold}.nc startfi.nc
else
   if (( ${numold} == 99999 )) ; then
    echo "No run because previous run crashed ! (99999 in num_run)"
    exit
   else
   echo "Where is file start"${numold}".nc??"
   exit
   fi
fi

# Run GCM -- THIS LINE NEEDS TO BE MODIFIED WITH THE CORRECT GCM EXECUTION COMMAND
mpirun -np 8 gcm_64x48x26_phystd_para.e < diagfi.def > lrun${numnew}


# Check if run ended normaly and copy datafiles
if [[ ( -r restartfi.nc  &&  -r restart.nc ) ]] ; then
  echo "Run seems to have ended normaly"


  \mv -f restart.nc start${numnew}.nc
  \mv -f restartfi.nc startfi${numnew}.nc  
    
else
  if [[ -r num_run ]] ; then
    \mv -f num_run num_run.crash
  else
    echo "No file num_run to build num_run.crash from !!"
    # Impose a default value of 0 for num_run
    echo 0 > num_run.crash
  fi
 echo 99999 > num_run
 exit
fi

# Copy other datafiles that may have been generated
if [[ -r diagfi.nc ]] ; then
  \mv -f diagfi.nc diagfi${numnew}.nc
fi
if [[ -r diagsoil.nc ]] ; then
  \mv -f diagsoil.nc diagsoil${numnew}.nc
fi
if [[ -r stats.nc ]] ; then
  \mv -f stats.nc stats${numnew}.nc
fi
if [[ -f profiles.dat ]] ; then
  \mv -f profiles.dat profiles${numnew}.dat
  \mv -f profiles.hdr profiles${numnew}.hdr
fi

# Prepare things for upcoming runs by writing
# value of computed season in file num_run
echo ${numnew} > num_run

# If we are over nummax : stop
if (( $numnew + 1 > $nummax )) ; then
   exit
# If not : restart the loop (copy the executable and run the copy)
else
   \cp -f run_gnome exe_mars
   ./exe_mars
fi

5

6

100

ncra -F -d Time,1,,1 diagfi.nc diagfi_MEAN.nc # format is "-d dimension,minimum,maximum,stride"

ncea -F -d Time,first,last diagfi.nc diagfi_subset.nc # format is "-d dimension,minimum,maximum" ; we recall you can type "ncdump -v time diagfi.nc" to see the Time values in the netCDF file.

ncks -F -d Time,1,,8 diagfi.nc diagfi_decimated.nc # format is "-d dimension,minimum,maximum,stride" ; In this example, this means that data is extracted 1 time every 8 time steps, starting from the first time step (number 1), ending at the last time step).

ncks -v tsurf,temp,p diagfi.nc diagfi_out.nc # Here we created a new file named 'diagfi_out.nc' in which we only kept variables named 'tsurf' (surface temperatures), 'temp' (atmospheric temperatures) and p (atmospheric pressures).

sudo apt install java
sudo apt install default-jre

panoply.sh

sudo apt install ncview

ncview diagfi.nc

 1 import numpy
 2 from netCDF4 import Dataset
 3 
 4 # HERE WE OPEN THE NETCDF FILE
 5 nc = Dataset('diagfi.nc')
 6 
 7 # HERE WE READ THE VARIABLES (1D OUTPUT)
 8 Time=nc.variables['Time'][:]
 9 lat=nc.variables['latitude'][:]
10 lon=nc.variables['longitude'][:]
11 al=nc.variables['altitude'][:]
12 
13 # HERE WE READ THE AREA (2D OUTPUT)
14 aire_GCM=nc.variables['aire'][:][:]
15 
16 # HERE WE READ 3D OUTPUTS
17 tsurf=nc.variables['tsurf'][:][:][:] # this is the surface temperature 3D field (time, latitude, longitude, altitude)
18 
19 # HERE WE READ 4D OUTPUTS
20 temp=nc.variables['temp'][:][:][:][:] # this is the atmospheric temperature 4D field (time, latitude, longitude, altitude)

1 from numpy import *
2 import numpy as np
3 
4 mean_tsurf=np.zeros((len(lat),len(lon)),dtype='f')
5 
6 for i in range(0,len(Time)):
7     for j in range(0,len(lat)):
8         for k in range(0,len(lon)):
9             mean_tsurf[j,k]=mean_tsurf[j,k]+tsurf[i,j,k]*(1./len(Time))

1 import matplotlib.pyplot as plt
2 
3 plt.figure(1)
4 plt.contourf(lon_GCM,lat_GCM,mean_tsurf)
5 plt.colorbar(label='Surface Temperature (K)')
6 plt.xlabel('Longitude ($^{\circ}$)')
7 plt.ylabel('Latitude ($^{\circ}$)')
8 plt.show()

 1 import numpy as np
 2 import xarray as xr 
 3 import matplotlib.pyplot as plt 
 4 
 5 # HERE WE OPEN THE NETCDF FILE
 6 data = xr.open_dataset('diagfi.nc',
 7                        decode_times=False)
 8 
 9 # HERE WE READ THE VARIABLES (1D OUTPUT)
10 Time=data['Time']
11 lat=data['latitude']
12 lon=data['longitude']
13 al=data['altitude']
14 
15 # HERE WE READ THE AREA (2D OUTPUT)
16 aire_GCM=data['aire']
17 
18 # HERE WE READ 3D OUTPUTS
19 tsurf=data['tsurf'] # this is the surface temperature 3D field (time, latitude, longitude, altitude)
20 
21 # HERE WE READ 4D OUTPUTS
22 temp=data['temp'] # this is the atmospheric temperature 4D field (time, latitude, longitude, altitude)
23 
24 ## let's take the time-averaged surface temperature
25 mean_tsurf = np.mean(tsurf,axis=0)
26 
27 ##Let's plot a lon-lat map
28 fig = plt.figure()
29 plt.contourf(lon,lat,mean_tsurf)
30 plt.colorbar(label='Surface Temperature (K)')
31 plt.xlabel('Longitude ($^{\circ}$)')
32 plt.ylabel('Latitude ($^{\circ}$)')
33 plt.show()

 1 import xarray as xr
 2 import os
 3 
 4 # your folder where output files are stored
 5 FOLDER = './your_folder_with_output_files/'
 6 
 7 # take back the files from your FOLDER
 8 list_files_folder=os.listdir(FOLDER)
 9 
10 # If there are several files.
11 # Sort your simulation files by date,
12 # so beginning of simulation will be top of the list
13 # and end of simulation will be end of the list.
14 list_files_folder.sort()
15 
16 files = [FOLDER+str(f) for f in list_files_folder]
17 # if you want to keep only files of special_year you can add this option :
18 # files = [FOLDER+str(f) for f in list_files_folder if f.startswith("special_year")]
19 
20 # xarray will magically concatenate your outfiles by 'Time' (or any other 'concat_dime' you want)
21 nc=xr.open_mfdataset(files,decode_times=False, concat_dim='Time', combine='nested')
22 
23 # to check your keys
24 for key in nc.keys():
25     print(key)
26 
27 # to load keys (example here with keys for a mesoscale simulation)
28 Times = nc['Times'][:]
29 PTOT = nc['PTOT'][:]
30 T = nc['T'][:]
31 W = nc['W'][:]
32 
33 # you can use some functions to make averages etc
34 
35 T_moy = T.mean(dim=['Time','south_north','west_east'])
36 
37 # other functions
38 # .cumsum (cumulative sum)
39 # .rename (change the name of the object)

 1 import xarray as xr
 2 import matplotlib.pyplot as plt 
 3 ##Load your data and print it 
 4 file = '/home/lteinturier/Documents/PhD/wasp43b/chemistry_project/input_5xsolar.nc'
 5 data = xr.open_dataset(file,decode_times=False)
 6 print(data)
 7 ##extracting the altitude level #20 for the whole file
 8 data = data.isel(altitude=20)
 9 ##let's assume that data hold a time-series of the temperature. Let's average it in time
10 temp = data['temp'].mean("Time",keep_attrs=True) #we keep the attribute when averaging to conserve the DataArray structure
11 ##now we plot 
12 fig = temp.plot.contourf(cmap='gnuplot',levels=50) #choose the colormap and the number of contourf levels 
13 fig.ax.set_title("P = {:.2e} mbar".format(data.p.mean().values)) ##set-up your title. If you don't change it, the title will be the altitude in km of your atmospheric level
14 plt.show()