#!/bin/env cdat # -*- coding: iso-8859-15 -*- # python system import os, sys, gc, getopt, time #import string import re import Numeric # CDAT/CDMS import cdms, MV #, cdutil from cdms.axis import TransientAxis from cdms.coord import TransientVirtualAxis from dataproc.utils import cdat def usage(): """Display usage information.""" man = """ This script with corrected the mask of a restart. It will get file_to_modify datas and the mask of file_to_add file and make a new file_out. The missing datas in file_to_modify will be simply replaced by the same points in file_to_add file. Usage : prompt> cdat ajoute_bords.py file_to_add file_to_modify file_out [variable] Add some points defines in file_to_add for all variables (or just one) of file_to_modify to file_out. examples : > ajoute_bords.py CM1414A_20001231_sechiba_rest.nc regrid_CF_CM5PIRC11_29901231_sechiba_rest.nc grille_regrid_CF_CM5PIRC11_29901231_sechiba_rest.nc > ajoute_bords.py CM1414A_20001231_stomate_rest.nc regrid_CF_CM5PIRC11_29901231_stomate_rest.nc grille_regrid_CF_CM5PIRC11_29901231_stomate_rest.nc > ajoute_bords.py CM5PIRC8_26601231_stomate_rest.nc regrid_CF_OOL-CTRLPI-1_9_5_SPIN_v2ORC_43_01001231_stomate_rest.nc grille_regrid_CF_OOL-CTRLPI-1_9_5_SPIN_v2ORC_43_01001231_stomate_rest.nc """ print man options = [ "help", "version" ] # Extract arguments from the command line try: myopts, myargs = getopt.getopt(sys.argv[1:], 'Vh', longopts = options) except getopt.GetoptError: print 'getopt error' usage() sys.exit(1) if (len(myargs) < 3): print 'I need at least three files.' usage() sys.exit(1) else: FileIn = myargs[0] print "file_to_add : ",FileIn if (os.path.exists(FileIn)==False): print 'file ',FileIn,' doesn''t exist !' sys.exit(1) else: fin=cdms.open(FileIn, 'r') FileOut = myargs[1] print "file_to_modify : ",FileOut if (os.path.exists(FileOut)==False): print 'file ',FileOut,' doesn''t exist !' sys.exit(1) else: fout=cdms.open(FileOut, 'r') FileNew = myargs[2] print "file_out : ",FileNew fnew=cdms.open(FileNew, 'w') if (len(myargs) > 3): varlist = [myargs[3]] print 'Only one variable :',varlist else: varlist = fin.listvariables() #varlist = [varlist[0]] print "Liste des variables de "+FileIn,varlist # Affichage des variables : fixe la limite max MV.set_print_limit(1700) var0=fin(varlist[0]) missing = var0.missing_value missing_ = var0.missing_value - 1. print missing, missing_ for var in varlist: print "-----" print "Variable : ",var,"------" if (var in ["nav_lon","nav_lat","nav_lev","time_steps","time", "lat","lon","lev","day_counter", "dt_days", "routingcounter", "veget_year", "date" ]): var1=fin(var) fnew.write(var1,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id=var1.id) print "remplace." continue if (var in ["bounds_lat","bounds_lon","routingcounter","routingarea","routetogrid","routetobasin", "basinid","topoindex","fastres","slowres","streamres","lakeres","previous_outflow", "lakeinflow","returnflow","riverflow","coastalflow","hydrographs","runoff_route", "drainage_route","evapot_route","precip_route","humrel_route","totnobio_route","vegtot_route"]): print "supprime." continue var0 = fin(var) var0m = MV.getmaskarray(var0) #var0_notmissing=(var0 < missing_) var0_notmissing=(var0m == 0) var1 = fout(var) var1m = MV.getmaskarray(var1) #var1_missing=(var1 > missing_) var1_missing=(var1m == 1) # var1or = MV.logical_or(var1,var1m) #test = MV.logical_and(var0_notmissing,var1or) test = MV.logical_and(var0_notmissing,var1_missing) # if ( var in ["temp_sol"] ): # fnew.write(var0m,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id="var0m") # fnew.write(var0_notmissing,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id="var0_notmissing") # fnew.write(var1_missing,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id="var1_missing") # fnew.write(var1m,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id="var1m") # fnew.write(var1or,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id="var1or") # fnew.write(test,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id="test") var1b = MV.where( test , var0, var1) # var1w = cdms.createVariable(var1b, # fill_value = missing, # attributes = cdat.getAttributes(var1), # id = var, # axes = var1.getAxisList()) # print var1w # var1w = cdms.createVariable(var1b, # attributes = cdat.getAttributes(var1), # id = var, # axes = [var1.getAxis(0), var1.getAxis(1), var1.getAxis(2)]) # print var1w,cdat.getAttributes(var1w), fnew.write(var1b,attributes=cdat.getAttributes(var1), axes=var1.getAxisList(), id=var1.id) fnew.close()