Changeset 5939


Ignore:
Timestamp:
10/11/21 14:30:35 (2 years ago)
Author:
snguyen
Message:

final version of Build_coordinates_mask.ipynb, tested with ORCA2.3, eORCA1.2 and PALEORCA

File:
1 edited

Legend:

Unmodified
Added
Removed
  • TOOLS/MOSAIX/Build_coordinates_mask.ipynb

    r5930 r5939  
    3333  { 
    3434   "cell_type": "markdown", 
     35   "id": "fb78f0e5", 
     36   "metadata": {}, 
     37   "source": [ 
     38    "# Various Switches To Deal With Bugs and Implementation Choices" 
     39   ] 
     40  }, 
     41  { 
     42   "cell_type": "code", 
     43   "execution_count": 2, 
     44   "id": "2bd430d4", 
     45   "metadata": {}, 
     46   "outputs": [], 
     47   "source": [ 
     48    "# reproduce periodicity bug for U and nperio = 4\n", 
     49    "nemoUbug = False\n", 
     50    "# change metrics (and bathymetry) for straits\n", 
     51    "straits = False\n", 
     52    "# periodicity imposed on coordinates, metrics and areas\n", 
     53    "coordperio = False\n", 
     54    "# function used to compute mask_T from bathymetry\n", 
     55    "maskbathynemo = False" 
     56   ] 
     57  }, 
     58  { 
     59   "cell_type": "markdown", 
    3560   "id": "26422e3f", 
    3661   "metadata": {}, 
     
    4166  { 
    4267   "cell_type": "code", 
    43    "execution_count": 2, 
     68   "execution_count": 3, 
    4469   "id": "deb66d3d", 
    4570   "metadata": {}, 
     
    4772   "source": [ 
    4873    "#model   = 'eORCA1.2'\n", 
    49     "model   = 'ORCA2.3'\n", 
    50     "#model   = 'ORCA2.4'" 
    51    ] 
    52   }, 
    53   { 
    54    "cell_type": "code", 
    55    "execution_count": 5, 
     74    "#model   = 'ORCA2.3'\n", 
     75    "model   = 'paleORCA'\n", 
     76    "#model   = 'ORCA2.4'\n", 
     77    "#model   = 'ORCA.05'http://localhost:8888/notebooks/Build_coordinates_mask_new_2409.ipynb#" 
     78   ] 
     79  }, 
     80  { 
     81   "cell_type": "code", 
     82   "execution_count": 4, 
    5683   "id": "2475c8b7", 
    5784   "metadata": {}, 
     
    5986   "source": [ 
    6087    "if model == 'eORCA1.2' :\n", 
    61     "    n_coord = 'eORCA_R1_coordinates_v1.0.nc' \n", 
    62     "    n_bathy = 'eORCA_R1_bathy_meter_v2.2.nc'\n", 
     88    "    n_coord = 'eORCA1.2coordinates.nc' # 'eORCA_R1_coordinates_v1.0.nc' \n", 
     89    "    n_bathy = 'eORCA1.2bathy_meter.nc' # 'eORCA_R1_bathy_meter_v2.2.nc'\n", 
    6390    "    #n_coord = 'https://vesg.ipsl.upmc.fr/thredds/dodsC/work/p86mart/GRAF/DATA/eORCA_R1_coordinates_v1.0.nc'\n", 
    6491    "    #n_bathy = 'https://vesg.ipsl.upmc.fr/thredds/dodsC/work/p86mart/GRAF/DATA/eORCA_R1_bathy_meter_v2.2.nc'\n", 
    6592    "    nperio  = 6\n", 
    66     "    n_out   = 'eORCA_R1_coordinates_mask_test.nc'\n", 
     93    "    n_out   = 'My_eORCA_R1_coordinates_mask_test.nc'\n", 
    6794    "\n", 
    6895    "if model == 'ORCA2.3' :\n", 
    6996    "    nperio  = 4\n", 
    70     "    n_coord = 'coordinates_ORCA2_LIM3_PISCES.nc' #'ORCA2.3_coordinates.nc' #'coordinates.nc' # \n", 
    71     "    n_bathy = 'bathy_meter_ORCA2_LIM3_PISCES.nc' #'ORCA2.3_bathy_meter.nc' #'bathy_meter_closea_topo.nc' # \n", 
     97    "    n_coord = 'coordinates.nc' # 'ORCA2.3_coordinates.nc' # 'coordinates_ORCA2_LIM3_PISCES.nc' #\n", 
     98    "    n_bathy = 'bathy_meter_closea_topo.nc' # 'ORCA2.3_bathy_meter.nc' # 'bathy_meter_ORCA2_LIM3_PISCES.nc' #\n", 
    7299    "    n_out   = 'My_test.nc'\n", 
    73100    "    \n", 
     
    76103    "    n_coord = 'ORCA2.4_coordinates.nc'\n", 
    77104    "    n_bathy = 'ORCA2.4_bathy_meter.nc'\n", 
    78     "    n_out   = 'ORCA2.4_coordinates_mask_test.nc'" 
     105    "    n_out   = 'ORCA2.4_coordinates_mask_test.nc'\n", 
     106    "    \n", 
     107    "if model == 'paleORCA' :\n", 
     108    "    n_coord = 'coordinates_paleo.nc'\n", 
     109    "    n_bathy = 'bathy_meter_paleo.nc'\n", 
     110    "    nperio = 6\n", 
     111    "    n_out = 'My_coordinates_mask_paleo.nc'\n", 
     112    "\n" 
    79113   ] 
    80114  }, 
     
    85119   "source": [ 
    86120    "## Open input files\n", 
    87     "Do not decode time which might be buggy in soem files\n", 
     121    "Do not decode time which might be buggy in some files\n", 
    88122    "Suppress singleton dimensions if necessary" 
    89123   ] 
     
    91125  { 
    92126   "cell_type": "code", 
    93    "execution_count": 6, 
     127   "execution_count": 5, 
    94128   "id": "0b7c7279", 
    95129   "metadata": { 
     
    104138  { 
    105139   "cell_type": "code", 
    106    "execution_count": 7, 
     140   "execution_count": 6, 
    107141   "id": "9fb92919", 
    108142   "metadata": {}, 
     
    112146     "output_type": "stream", 
    113147     "text": [ 
    114       "success\n" 
     148      "failed to suppress time\n" 
    115149     ] 
    116150    } 
     
    144178  { 
    145179   "cell_type": "code", 
    146    "execution_count": 8, 
     180   "execution_count": 7, 
    147181   "id": "9dfd90e6", 
    148182   "metadata": { 
     
    154188     "output_type": "stream", 
    155189     "text": [ 
    156       "Exception\n" 
     190      "Normal case\n" 
    157191     ] 
    158192    } 
    159193   ], 
    160194   "source": [ 
    161     "Bathymetry = f_bathy['Bathymetry']\n", 
     195    "Bathymetry = f_bathy['Bathymetry'].copy()\n", 
    162196    "\n", 
    163197    "try :\n", 
     
    176210  { 
    177211   "cell_type": "code", 
     212   "execution_count": 8, 
     213   "id": "2606d14b", 
     214   "metadata": {}, 
     215   "outputs": [], 
     216   "source": [ 
     217    "if straits :\n", 
     218    "    # Open straits for ORCA2 grid\n", 
     219    "    if model == 'ORCA2.3' :\n", 
     220    "        # orca_r2: Gibraltar strait open\n", 
     221    "        Bathymetry[101,139] = 284.\n", 
     222    "        # orca_r2: Bab el Mandeb strait open\n", 
     223    "        Bathymetry[87,159] = 137." 
     224   ] 
     225  }, 
     226  { 
     227   "cell_type": "markdown", 
     228   "id": "71c9c821", 
     229   "metadata": {}, 
     230   "source": [ 
     231    "### Determine which periodicity is needed" 
     232   ] 
     233  }, 
     234  { 
     235   "cell_type": "code", 
    178236   "execution_count": 9, 
    179    "id": "2606d14b", 
    180    "metadata": {}, 
    181    "outputs": [], 
    182    "source": [ 
    183     "if model == 'ORCA2.3' :\n", 
    184     "    # orca_r2: Gibraltar strait open\n", 
    185     "    Bathymetry[101,139] = 284.\n", 
    186     "    # orca_r2: Bab el Mandeb strait open\n", 
    187     "    Bathymetry[87,159] = 137." 
    188    ] 
    189   }, 
    190   { 
    191    "cell_type": "markdown", 
    192    "id": "71c9c821", 
    193    "metadata": {}, 
    194    "source": [ 
    195     "### Determine which periodicity is needed" 
    196    ] 
    197   }, 
    198   { 
    199    "cell_type": "code", 
    200    "execution_count": 10, 
    201237   "id": "7e7a49d5", 
    202238   "metadata": {}, 
     
    206242     "output_type": "stream", 
    207243     "text": [ 
    208       "nperio specified :  4\n" 
     244      "nperio specified :  6\n" 
    209245     ] 
    210246    } 
     
    229265  { 
    230266   "cell_type": "code", 
     267   "execution_count": 10, 
     268   "id": "34cdcfe4", 
     269   "metadata": {}, 
     270   "outputs": [], 
     271   "source": [ 
     272    "Bathymetry = nemo.lbc (Bathymetry, nperio=nperio, cd_type='T')\n", 
     273    "#Suppress ocean points at southernmost position\n", 
     274    "Bathymetry[0,:] = 0.0 " 
     275   ] 
     276  }, 
     277  { 
     278   "cell_type": "markdown", 
     279   "id": "0da91538", 
     280   "metadata": {}, 
     281   "source": [ 
     282    "### Creates ```mask_T```" 
     283   ] 
     284  }, 
     285  { 
     286   "cell_type": "code", 
    231287   "execution_count": 11, 
    232    "id": "34cdcfe4", 
    233    "metadata": {}, 
    234    "outputs": [], 
    235    "source": [ 
    236     "Bathymetry = nemo.lbc (Bathymetry, nperio=nperio, cd_type='T')" 
    237    ] 
    238   }, 
    239   { 
    240    "cell_type": "markdown", 
    241    "id": "0da91538", 
    242    "metadata": {}, 
    243    "source": [ 
    244     "### Creates ```mask_T```" 
     288   "id": "49a654ef", 
     289   "metadata": {}, 
     290   "outputs": [], 
     291   "source": [ 
     292    "if maskbathynemo :\n", 
     293    "    # Use same formula as domzgr. Should be equivalent to Bathimetry > 0.0. \n", 
     294    "    mask_T = xr.where (Bathymetry - 1. + 0.1 >= 0.0, 1, 0).astype (dtype='f4')\n", 
     295    "else :\n", 
     296    "    mask_T = xr.where (Bathymetry > 0.0, 1, 0).astype (dtype='f4')" 
     297   ] 
     298  }, 
     299  { 
     300   "cell_type": "markdown", 
     301   "id": "327d48ea", 
     302   "metadata": {}, 
     303   "source": [ 
     304    "### Creates other masks" 
    245305   ] 
    246306  }, 
     
    248308   "cell_type": "code", 
    249309   "execution_count": 12, 
    250    "id": "49a654ef", 
    251    "metadata": {}, 
    252    "outputs": [], 
    253    "source": [ 
    254     "mask_T = xr.where (Bathymetry - 1. + 0.1 >= 0.0, 1, 0).astype (dtype='int8')" 
    255    ] 
    256   }, 
    257   { 
    258    "cell_type": "markdown", 
    259    "id": "327d48ea", 
    260    "metadata": {}, 
    261    "source": [ 
    262     "### Creates other masks" 
    263    ] 
    264   }, 
    265   { 
    266    "cell_type": "code", 
    267    "execution_count": 14, 
    268310   "id": "f0b656a1", 
    269311   "metadata": {}, 
     
    275317    "mask_W = mask_T\n", 
    276318    "\n", 
    277     "mask_U = nemo.lbc (mask_U, nperio=nperio, cd_type='U', nemo_4U_bug=True).astype (dtype='int8')\n", 
    278     "mask_V = nemo.lbc (mask_V, nperio=nperio, cd_type='V').astype (dtype='int8')\n", 
    279     "mask_F = nemo.lbc (mask_F, nperio=nperio, cd_type='F').astype (dtype='int8')\n", 
    280     "mask_W = nemo.lbc (mask_W, nperio=nperio, cd_type='W').astype (dtype='int8')\n", 
     319    "mask_U = nemo.lbc (mask_U, nperio=nperio, cd_type='U', nemo_4U_bug=nemoUbug).astype (dtype='f4')\n", 
     320    "mask_V = nemo.lbc (mask_V, nperio=nperio, cd_type='V').astype (dtype='f4')\n", 
     321    "mask_F = nemo.lbc (mask_F, nperio=nperio, cd_type='F').astype (dtype='f4')\n", 
     322    "mask_W = nemo.lbc (mask_W, nperio=nperio, cd_type='W').astype (dtype='f4')\n", 
    281323    "\n", 
    282324    "mask_U = mask_U.rename ( {'y_grid_T'      : 'y_grid_U'      , 'x_grid_T'      : 'x_grid_U', \n", 
     
    306348  { 
    307349   "cell_type": "code", 
    308    "execution_count": null, 
    309    "id": "b9ce5946", 
    310    "metadata": {}, 
    311    "outputs": [], 
    312    "source": [ 
    313     "#maskutil_T = nemo.lbc_mask (mask_T, nperio=nperio, cd_type='T')\n", 
    314     "#maskutil_U = nemo.lbc_mask (mask_U, nperio=nperio, cd_type='U')\n", 
    315     "#maskutil_V = nemo.lbc_mask (mask_V, nperio=nperio, cd_type='V')\n", 
    316     "#maskutil_F = nemo.lbc_mask (mask_F, nperio=nperio, cd_type='F')\n", 
    317     "#maskutil_W = nemo.lbc_mask (mask_W, nperio=nperio, cd_type='W')\n", 
    318     "\n", 
    319     "#maskutil_T.name = 'maskutil_T'\n", 
    320     "#maskutil_U.name = 'maskutil_U'\n", 
    321     "#maskutil_V.name = 'maskutil_V'\n", 
    322     "#maskutil_F.name = 'maskutil_F'\n", 
    323     "#maskutil_W.name = 'maskutil_W'" 
    324    ] 
    325   }, 
    326   { 
    327    "cell_type": "code", 
    328    "execution_count": 15, 
     350   "execution_count": 13, 
    329351   "id": "274e4321", 
    330352   "metadata": {}, 
    331353   "outputs": [], 
    332354   "source": [ 
    333     "# Plus compact que ci-dessus, mais un peu incompréhensible .... !!!\n", 
    334355    "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 
    335356    "    MaskName = 'mask_'     + cd_type\n", 
    336357    "    UtilName = 'maskutil_' + cd_type\n", 
    337     "    locals()[UtilName] = nemo.lbc_mask (locals()[MaskName], nperio=nperio, cd_type=cd_type)\n", 
    338     "    locals()[UtilName].name = UtilName" 
     358    "    locals()[UtilName] = nemo.lbc_mask (locals()[MaskName].copy(), nperio=nperio, cd_type=cd_type).astype (dtype='f4')\n", 
     359    "    locals()[UtilName].name = UtilName\n", 
     360    "    locals()[UtilName].encoding['_FillValue'] = None" 
    339361   ] 
    340362  }, 
     
    349371  { 
    350372   "cell_type": "code", 
    351    "execution_count": 16, 
     373   "execution_count": 14, 
    352374   "id": "f46e1c93", 
    353375   "metadata": {}, 
     
    377399  { 
    378400   "cell_type": "code", 
    379    "execution_count": 17, 
     401   "execution_count": 15, 
    380402   "id": "a7ce1490", 
    381403   "metadata": {}, 
    382404   "outputs": [], 
    383405   "source": [ 
    384     "nav_lon_grid_T = nemo.lbc (f_coord ['glamt'], nperio=nperio, cd_type='T')\n", 
    385     "nav_lat_grid_T = nemo.lbc (f_coord ['gphit'], nperio=nperio, cd_type='T')\n", 
    386     "nav_lon_grid_U = nemo.lbc (f_coord ['glamu'], nperio=nperio, cd_type='U',nemo_4U_bug=True)\n", 
    387     "#nav_lon_grid_U = f_coord ['glamu']\n", 
    388     "nav_lat_grid_U = nemo.lbc (f_coord ['gphiu'], nperio=nperio, cd_type='U',nemo_4U_bug=True)\n", 
    389     "#nav_lat_grid_U = f_coord ['gphiu']\n", 
    390     "nav_lon_grid_V = nemo.lbc (f_coord ['glamv'], nperio=nperio, cd_type='V')\n", 
    391     "nav_lat_grid_V = nemo.lbc (f_coord ['gphiv'], nperio=nperio, cd_type='V')\n", 
    392     "nav_lon_grid_W = nemo.lbc (f_coord ['glamt'], nperio=nperio, cd_type='W')\n", 
    393     "nav_lat_grid_W = nemo.lbc (f_coord ['gphit'], nperio=nperio, cd_type='W')\n", 
    394     "nav_lon_grid_F = nemo.lbc (f_coord ['glamf'], nperio=nperio, cd_type='F')\n", 
    395     "nav_lat_grid_F = nemo.lbc (f_coord ['gphif'], nperio=nperio, cd_type='F')\n", 
     406    "if coordperio :\n", 
     407    "    nav_lon_grid_T = nemo.lbc (f_coord ['glamt'].copy(), nperio=nperio, cd_type='T')\n", 
     408    "    nav_lat_grid_T = nemo.lbc (f_coord ['gphit'].copy(), nperio=nperio, cd_type='T')\n", 
     409    "    nav_lon_grid_U = nemo.lbc (f_coord ['glamu'].copy(), nperio=nperio, cd_type='U',nemo_4U_bug=nemoUbug)\n", 
     410    "    nav_lat_grid_U = nemo.lbc (f_coord ['gphiu'].copy(), nperio=nperio, cd_type='U',nemo_4U_bug=nemoUbug)\n", 
     411    "    nav_lon_grid_V = nemo.lbc (f_coord ['glamv'].copy(), nperio=nperio, cd_type='V')\n", 
     412    "    nav_lat_grid_V = nemo.lbc (f_coord ['gphiv'].copy(), nperio=nperio, cd_type='V')\n", 
     413    "    nav_lon_grid_W = nemo.lbc (f_coord ['glamt'].copy(), nperio=nperio, cd_type='W')\n", 
     414    "    nav_lat_grid_W = nemo.lbc (f_coord ['gphit'].copy(), nperio=nperio, cd_type='W')\n", 
     415    "    nav_lon_grid_F = nemo.lbc (f_coord ['glamf'].copy(), nperio=nperio, cd_type='F')\n", 
     416    "    nav_lat_grid_F = nemo.lbc (f_coord ['gphif'].copy(), nperio=nperio, cd_type='F')\n", 
     417    "else :\n", 
     418    "    nav_lon_grid_T = f_coord ['glamt'].copy()\n", 
     419    "    nav_lat_grid_T = f_coord ['gphit'].copy()\n", 
     420    "    nav_lon_grid_U = f_coord ['glamu'].copy()\n", 
     421    "    nav_lat_grid_U = f_coord ['gphiu'].copy()\n", 
     422    "    nav_lon_grid_V = f_coord ['glamv'].copy()\n", 
     423    "    nav_lat_grid_V = f_coord ['gphiv'].copy()\n", 
     424    "    nav_lon_grid_W = f_coord ['glamt'].copy()\n", 
     425    "    nav_lat_grid_W = f_coord ['gphit'].copy()\n", 
     426    "    nav_lon_grid_F = f_coord ['glamf'].copy()\n", 
     427    "    nav_lat_grid_F = f_coord ['gphif'].copy()\n", 
    396428    "\n", 
    397429    "nav_lon_grid_T = nav_lon_grid_T.rename( {'y':'y_grid_T', 'x':'x_grid_T'} )\n", 
     
    406438    "nav_lat_grid_F = nav_lat_grid_F.rename( {'y':'y_grid_F', 'x':'x_grid_F'} )\n", 
    407439    "\n", 
     440    "# remove _FillValue and missing_value\n", 
    408441    "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 
    409442    "    for dir_type in ['lat', 'lon'] :\n", 
     
    415448  { 
    416449   "cell_type": "code", 
    417    "execution_count": 18, 
     450   "execution_count": 16, 
    418451   "id": "516b4c5c", 
    419452   "metadata": {}, 
     
    442475  { 
    443476   "cell_type": "code", 
    444    "execution_count": 27, 
     477   "execution_count": 17, 
    445478   "id": "de0e6514", 
    446479   "metadata": {}, 
     
    496529  { 
    497530   "cell_type": "code", 
    498    "execution_count": 20, 
     531   "execution_count": 18, 
    499532   "id": "8e4d997c", 
    500533   "metadata": {}, 
    501534   "outputs": [], 
    502535   "source": [ 
     536    "# create new variables e1 e2 to keep f_coord the same\n", 
     537    "for cd_type in ['t', 'u', 'v', 'f'] :\n", 
     538    "    for axis in ['1', '2'] :\n", 
     539    "        coordName = 'e' + axis + cd_type\n", 
     540    "        locals()[coordName]=f_coord[coordName].copy()\n", 
    503541    "# remove zero values from areas\n", 
    504542    "# need to be define for the extended grid south of -80S\n", 
    505543    "# some point are undefined but you need to have e1 and e2 .NE. 0\n", 
    506     "for cd_type in ['t', 'u', 'v', 'f'] :\n", 
    507     "    for axis in ['1', '2'] :\n", 
    508     "        coordName = 'e' + axis + cd_type\n", 
    509     "        f_coord[coordName]=xr.where(f_coord[coordName] == 0.0, 1.0e2, f_coord[coordName])" 
    510    ] 
    511   }, 
    512   { 
    513    "cell_type": "code", 
    514    "execution_count": 21, 
     544    "        locals()[coordName]=xr.where(locals()[coordName] == 0.0, 1.0e2, locals()[coordName])" 
     545   ] 
     546  }, 
     547  { 
     548   "cell_type": "code", 
     549   "execution_count": 19, 
    515550   "id": "175bae75", 
    516551   "metadata": {}, 
     
    518553   "source": [ 
    519554    "# Correct areas for straits\n", 
    520     "if model == 'ORCA2.3' :\n", 
    521     "    # orca_r2: Gibraltar    : e2u reduced to 20 km\n", 
    522     "    f_coord['e2u'][101,138:140] = 20.e3\n", 
    523     "    # orca_r2: Bab el Mandeb: e2u reduced to 30 km\n", 
    524     "    #                         e1v reduced to 18 km\n", 
    525     "    f_coord['e1v'][87,159] = 18.e3\n", 
    526     "    f_coord['e2u'][87,159] = 30.e3\n", 
    527     "    # orca_r2: Danish Straits : e2u reduced to 10 km\n", 
    528     "    f_coord['e2u'][115,144:146] = 10.e3\n", 
    529     "    " 
    530    ] 
    531   }, 
    532   { 
    533    "cell_type": "code", 
    534    "execution_count": 22, 
     555    "if straits :\n", 
     556    "    # ORCA R2 configuration\n", 
     557    "    if model == 'ORCA2.3' :\n", 
     558    "        # Gibraltar     : e2u reduced to 20 km\n", 
     559    "        e2u[101,138:140] = 20.e3\n", 
     560    "        # Bab el Mandeb : e2u reduced to 30 km\n", 
     561    "        #                e1v reduced to 18 km\n", 
     562    "        e1v[87,159] = 18.e3\n", 
     563    "        e2u[87,159] = 30.e3\n", 
     564    "        # Danish Straits: e2u reduced to 10 km\n", 
     565    "        e2u[115,144:146] = 10.e3\n", 
     566    "    # ORCA R1 configuration\n", 
     567    "    if model == 'eORCA1.2' :\n", 
     568    "        # Gibraltar : e2u reduced to 20 km\n", 
     569    "        e2u[240,281:283] = 20.e3\n", 
     570    "        # Bhosporus : e2u reduced to 10 km\n", 
     571    "        e2u[247,313:315] =  10.e3\n", 
     572    "        # Lombok : e1v reduced to 13 km\n", 
     573    "        e1v[163:165,43] =  13.e3\n", 
     574    "        # Sumba : e1v reduced to 8 km\n", 
     575    "        e1v[163:165,47] =  8.e3\n", 
     576    "        # Ombai : e1v reduced to 13 km\n", 
     577    "        e1v[163:165,52] = 13.e3\n", 
     578    "        # Timor Passage : e1v reduced to 20 km\n", 
     579    "        #e1v[163:165,55] = 20.e3\n", 
     580    "        # W Halmahera : e1v reduced to 30 km\n", 
     581    "        e1v[180:182,54] = 30.e3\n", 
     582    "        # E Halmahera : e1v reduced to 50 km\n", 
     583    "        e1v[180:182,57] = 50.e3\n", 
     584    "    # ORCA R05 configuration\n", 
     585    "    if model == 'ORCA.05' :\n", 
     586    "        # Reduced e2u at the Gibraltar Strait\n", 
     587    "        e2u[326,562:564] =  20.e3\n", 
     588    "        # Reduced e2u at the Bosphore Strait\n", 
     589    "        e2u[342,626:628] =  10.e3\n", 
     590    "        # Reduced e2u at the Sumba Strait\n", 
     591    "        e2u[231,92:94] =  40.e3\n", 
     592    "        # Reduced e2u at the Ombai Strait\n", 
     593    "        e2u[231,102] =  15.e3\n", 
     594    "        # Reduced e2u at the Palk Strait\n", 
     595    "        e2u[269,14] =  10.e3\n", 
     596    "        # Reduced e1v at the Lombok Strait\n", 
     597    "        e1v[231:233,86] =  10.e3\n", 
     598    "        # Reduced e1v at the Bab el Mandeb\n", 
     599    "        e1v[275,661] =  25.e3" 
     600   ] 
     601  }, 
     602  { 
     603   "cell_type": "code", 
     604   "execution_count": 20, 
    535605   "id": "4f6cfbca", 
    536606   "metadata": {}, 
    537607   "outputs": [], 
    538608   "source": [ 
    539     "area_grid_T = nemo.lbc (f_coord ['e1t']*f_coord ['e2t'], nperio=nperio, cd_type='T')\n", 
    540     "area_grid_U = nemo.lbc (f_coord ['e1u']*f_coord ['e2u'], nperio=nperio, cd_type='U',nemo_4U_bug=True)\n", 
    541     "#area_grid_U = f_coord ['e1u']*f_coord ['e2u']\n", 
    542     "area_grid_V = nemo.lbc (f_coord ['e1v']*f_coord ['e2v'], nperio=nperio, cd_type='V')\n", 
    543     "area_grid_W = nemo.lbc (f_coord ['e1t']*f_coord ['e2t'], nperio=nperio, cd_type='W')\n", 
    544     "area_grid_F = nemo.lbc (f_coord ['e1f']*f_coord ['e2f'], nperio=nperio, cd_type='F')\n", 
     609    "if coordperio :\n", 
     610    "    area_grid_T = nemo.lbc (e1t*e2t, nperio=nperio, cd_type='T')\n", 
     611    "    area_grid_U = nemo.lbc (e1u*e2u, nperio=nperio, cd_type='U',nemo_4U_bug=nemoUbug)\n", 
     612    "    area_grid_V = nemo.lbc (e1v*e2v, nperio=nperio, cd_type='V')\n", 
     613    "    area_grid_W = nemo.lbc (e1t*e2t, nperio=nperio, cd_type='W')\n", 
     614    "    area_grid_F = nemo.lbc (e1f*e2f, nperio=nperio, cd_type='F')\n", 
     615    "else :\n", 
     616    "    area_grid_T = e1t*e2t\n", 
     617    "    area_grid_U = e1u*e2u\n", 
     618    "    area_grid_V = e1v*e2v\n", 
     619    "    area_grid_W = e1t*e2t\n", 
     620    "    area_grid_F = e1f*e2f\n", 
     621    "\n", 
    545622    "\n", 
    546623    "area_grid_T.name = 'area_grid_T'\n", 
     
    558635    "area_grid_T.attrs['standard_name'] = 'cell_area'\n", 
    559636    "area_grid_T.attrs['units']         = 'm2'\n", 
    560     "area_grid_T.attrs['coordinates']   = 'nav_lat_grid_T nav_lon_grid_T'\n", 
     637    "#area_grid_T.attrs['coordinates']   = 'nav_lat_grid_T nav_lon_grid_T'\n", 
    561638    "area_grid_U.attrs['standard_name'] = 'cell_area'\n", 
    562639    "area_grid_U.attrs['units']         = 'm2'\n", 
    563     "area_grid_U.attrs['coordinates']   = 'nav_lat_grid_U nav_lon_grid_U'\n", 
     640    "#area_grid_U.attrs['coordinates']   = 'nav_lat_grid_U nav_lon_grid_U'\n", 
    564641    "area_grid_V.attrs['standard_name'] = 'cell_area'\n", 
    565642    "area_grid_V.attrs['units']         = 'm2'\n", 
    566     "area_grid_V.attrs['coordinates']   = 'nav_lat_grid_V nav_lon_grid_V'\n", 
     643    "#area_grid_V.attrs['coordinates']   = 'nav_lat_grid_V nav_lon_grid_V'\n", 
    567644    "area_grid_W.attrs['standard_name'] = 'cell_area'\n", 
    568645    "area_grid_W.attrs['units']         = 'm2'\n", 
    569     "area_grid_W.attrs['coordinates']   = 'nav_lat_grid_W nav_lon_grid_W'\n", 
     646    "#area_grid_W.attrs['coordinates']   = 'nav_lat_grid_W nav_lon_grid_W'\n", 
    570647    "area_grid_F.attrs['standard_name'] = 'cell_area'\n", 
    571648    "area_grid_F.attrs['units']         = 'm2'\n", 
    572     "area_grid_F.attrs['coordinates']   = 'nav_lat_grid_F nav_lon_grid_F'" 
     649    "#area_grid_F.attrs['coordinates']   = 'nav_lat_grid_F nav_lon_grid_F'\n", 
     650    "\n", 
     651    "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 
     652    "    areaName = 'area_grid_' + cd_type\n", 
     653    "    locals()[areaName].encoding['_FillValue'] = None" 
    573654   ] 
    574655  }, 
     
    583664  { 
    584665   "cell_type": "code", 
    585    "execution_count": 23, 
     666   "execution_count": 21, 
    586667   "id": "69c404d9", 
    587668   "metadata": {}, 
     
    597678    "    if cdgrd in ['T', 'W'] : \n", 
    598679    "        icnr = -1 ; jcnr = -1\n", 
    599     "        corner_lon = nav_lon_grid_F ; corner_lat = nav_lat_grid_F\n", 
     680    "        corner_lon = f_coord ['glamf'].copy() ; corner_lat = f_coord ['gphif'].copy()\n", 
     681    "        center_lon = f_coord ['glamt'].copy() ; center_lat = f_coord ['gphit'].copy()\n", 
    600682    "    if cdgrd == 'U'        : \n", 
    601683    "        icnr =  0 ; jcnr = -1\n", 
    602     "        corner_lon = nav_lon_grid_V ; corner_lat = nav_lat_grid_V\n", 
     684    "        corner_lon = f_coord ['glamv'].copy() ; corner_lat = f_coord ['gphiv'].copy()\n", 
     685    "        center_lon = f_coord ['glamu'].copy() ; center_lat = f_coord ['gphiu'].copy()\n", 
    603686    "    if cdgrd == 'V'        : \n", 
    604687    "        icnr = -1 ; jcnr =  0\n", 
    605     "        corner_lon = nav_lon_grid_U ; corner_lat = nav_lat_grid_U\n", 
     688    "        corner_lon = f_coord ['glamu'].copy() ; corner_lat = f_coord ['gphiu'].copy()\n", 
     689    "        center_lon = f_coord ['glamv'].copy() ; center_lat = f_coord ['gphiv'].copy()\n", 
    606690    "    if cdgrd == 'F'        : \n", 
    607691    "        icnr = -1 ; jcnr = -1\n", 
    608     "        corner_lon = nav_lon_grid_T ; corner_lat = nav_lat_grid_T\n", 
    609     "\n", 
    610     "    y_grid, x_grid = corner_lon.shape ; nvertex = 4\n", 
     692    "        corner_lon = f_coord ['glamt'].copy() ; corner_lat = f_coord ['gphit'].copy()\n", 
     693    "        center_lon = f_coord ['glamf'].copy() ; center_lat = f_coord ['gphif'].copy()\n", 
     694    "      \n", 
     695    "    #y_grid, x_grid = corner_lon.shape ;\n", 
     696    "    nvertex = 4\n", 
    611697    "    dims = ['y_grid_' + cdgrd, 'x_grid_' + cdgrd, 'nvertex_grid_' + cdgrd]\n", 
    612698    "    \n", 
    613     "    bounds_lon = xr.DataArray (np.zeros ((y_grid, x_grid, nvertex)), dims=dims)\n", 
    614     "    bounds_lat = xr.DataArray (np.zeros ((y_grid, x_grid, nvertex)), dims=dims)\n", 
     699    "    bounds_lon = xr.DataArray (np.zeros ((jpj, jpi, nvertex)), dims=dims)\n", 
     700    "    bounds_lat = xr.DataArray (np.zeros ((jpj, jpi, nvertex)), dims=dims)\n", 
    615701    "    \n", 
    616702    "    idx = [(jcnr,icnr), (jcnr,icnr+1), (jcnr+1,icnr+1), (jcnr+1,icnr)]\n", 
     
    619705    "    # and complete with periodicity\n", 
    620706    "    for nn in range (nvertex) :\n", 
    621     "        bounds_lon[:,:,nn] = np.roll  (corner_lon, shift=(idx[nn]), axis=(-2,-1))\n", 
    622     "        bounds_lat[:,:,nn] = np.roll  (corner_lat, shift=(idx[nn]), axis=(-2,-1))\n", 
    623     "        bounds_lon[:,:,nn] = nemo.lbc (bounds_lon[:,:,nn], nperio=nperio, cd_type=cdgrd)\n", 
    624     "        bounds_lat[:,:,nn] = nemo.lbc (bounds_lat[:,:,nn], nperio=nperio, cd_type=cdgrd)\n", 
     707    "        tmp = np.roll (corner_lon, shift=tuple(-1*np.array(idx[nn])), axis=(-2,-1))\n", 
     708    "        bounds_lon[1:jpj,1:jpi,nn] = tmp[1:jpj,1:jpi]\n", 
     709    "        tmp = np.roll (corner_lat, shift=tuple(-1*np.array(idx[nn])), axis=(-2,-1))\n", 
     710    "        bounds_lat[1:jpj,1:jpi,nn] = tmp[1:jpj,1:jpi]\n", 
     711    "        bounds_lon[:,:,nn] = nemo.lbc (bounds_lon[:,:,nn], nperio=nperio, cd_type=cdgrd, nemo_4U_bug=nemoUbug)\n", 
     712    "        bounds_lat[:,:,nn] = nemo.lbc (bounds_lat[:,:,nn], nperio=nperio, cd_type=cdgrd, nemo_4U_bug=nemoUbug)\n", 
     713    "    \n", 
     714    "    # Zero-size cells at closed boundaries if cell points provided,\n", 
     715    "    # otherwise they are closed cells with unrealistic bounds\n", 
     716    "    if not (nperio == 1 or nperio == 4 or nperio == 6) :\n", 
     717    "        for nn in range (nvertex) : \n", 
     718    "            bounds_lon[:,0,nn]   = center_lon[:,0]   # (West or jpni = 1), closed E-W\n", 
     719    "            bounds_lat[:,0,nn]   = center_lat[:,0]\n", 
     720    "            bounds_lon[:,jpi,nn] = center_lon[:,jpi] # (East or jpni = 1), closed E-W\n", 
     721    "            bounds_lat[:,jpi,nn] = center_lat[:,jpi]\n", 
     722    "    if nperio != 2 :\n", 
     723    "        for nn in range (nvertex) :\n", 
     724    "            bounds_lon[0,:,nn] = center_lon[0,:] # (South or jpnj = 1), not symmetric\n", 
     725    "            bounds_lat[0,:,nn] = center_lat[0,:]\n", 
     726    "    if nperio < 3 :\n", 
     727    "        for nn in range (nvertex) :\n", 
     728    "            bounds_lon[jpj,:,nn] = center_lon[jpj,:] # (North or jpnj = 1), no north fold\n", 
     729    "            bounds_lat[jpj,:,nn] = center_lat[jpj,:]\n", 
    625730    "    \n", 
    626731    "    # Rotate cells at the north fold\n", 
    627732    "    if nperio >= 3 :\n", 
    628733    "        # Working array for location of northfold\n", 
    629     "        z_fld = nemo.lbc (np.ones ((jpj, jpi)), nperio=nperio, cd_type=cdgrd, psgn=-1. )\n", 
    630     "        z_fld = (z_fld == -1.0)\n", 
    631     "        \n", 
    632     "        for nn in range (nvertex) :\n", 
    633     "            nr = np.mod (nn+2, 4)\n", 
    634     "            bounds_lon[:,:,nn] = np.where (z_fld, bounds_lon[:,:,nr], bounds_lon[:,:,nn])\n", 
    635     "            bounds_lat[:,:,nn] = np.where (z_fld, bounds_lat[:,:,nr], bounds_lat[:,:,nn])\n", 
     734    "        z_fld = nemo.lbc (np.ones ((jpj, jpi)), nperio=nperio, cd_type=cdgrd, psgn=-1., nemo_4U_bug=nemoUbug)\n", 
     735    "        z_fld = np.repeat((z_fld == -1.0)[...,np.newaxis],4,axis=2)\n", 
     736    "        # circular shift of 2 indices in bounds third index\n", 
     737    "        bounds_lon_tmp = np.roll (bounds_lon, shift=-2, axis=2)\n", 
     738    "        bounds_lat_tmp = np.roll (bounds_lat, shift=-2, axis=2)\n", 
     739    "        bounds_lon[:,:,:]  = np.where (z_fld, bounds_lon_tmp[:,:,:] , bounds_lon[:,:,:] )\n", 
     740    "        bounds_lat[:,:,:]  = np.where (z_fld, bounds_lat_tmp[:,:,:] , bounds_lat[:,:,:] )\n", 
    636741    "    \n", 
    637742    "    # Invert cells at the symmetric equator\n", 
    638743    "    if nperio == 2 :\n", 
    639     "        for nn in range (nvertex) :\n", 
    640     "            nr = np.mod (nn+2, 4)\n", 
    641     "            bounds_lon [0, :, nn] = bounds_lon [0, :, nr]\n", 
    642     "            bounds_lat [0, :, nn] = bounds_lat [0, :, nr]\n", 
    643     "    \n", 
    644     "    bounds_lon.attrs['coordinates'] = 'nav_lat_grid_' + cdgrd + ' nav_lon_grid_' + cdgrd\n", 
    645     "    bounds_lat.attrs['coordinates'] = 'nav_lat_grid_' + cdgrd + ' nav_lon_grid_' + cdgrd\n", 
     744    "        bounds_lon_tmp = np.roll (bounds_lon, shift=-2, axis=2)\n", 
     745    "        bounds_lat_tmp = np.roll (bounds_lat, shift=-2, axis=2)\n", 
     746    "        bounds_lon [0,:,:] = bounds_lon [0,:,:]\n", 
     747    "        bounds_lat [0,:,:] = bounds_lat [0,:,:]\n", 
     748    "    \n", 
     749    "    #bounds_lon.attrs['coordinates'] = 'nav_lat_grid_' + cdgrd + ' nav_lon_grid_' + cdgrd\n", 
     750    "    #bounds_lat.attrs['coordinates'] = 'nav_lat_grid_' + cdgrd + ' nav_lon_grid_' + cdgrd\n", 
    646751    "    bounds_lon.attrs['units']       = 'degrees_east'\n", 
    647752    "    bounds_lat.attrs['units']       = 'degrees_north'\n", 
    648753    "    bounds_lon.name = 'bounds_lon_grid_' + cdgrd\n", 
    649754    "    bounds_lat.name = 'bounds_lat_grid_' + cdgrd\n", 
     755    "    # remove _FillValue\n", 
     756    "    bounds_lon.encoding['_FillValue'] = None\n", 
     757    "    bounds_lat.encoding['_FillValue'] = None\n", 
    650758    "\n", 
    651759    "    return bounds_lon, bounds_lat" 
     
    654762  { 
    655763   "cell_type": "code", 
    656    "execution_count": 24, 
     764   "execution_count": 22, 
    657765   "id": "87a8f592", 
    658766   "metadata": {}, 
     
    676784  { 
    677785   "cell_type": "code", 
    678    "execution_count": 28, 
     786   "execution_count": 23, 
    679787   "id": "ca87fa1b", 
    680788   "metadata": {}, 
     
    708816    "    'bounds_lat_grid_F': bounds_lat_grid_F,\n", 
    709817    "})\n", 
     818    "#replace nav_lon nav_lat with variables obtained from NEMO coordinates.nc file\n", 
     819    "#by construction nav_lon nav_lat are those of the bathymetry \n", 
    710820    "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 
    711821    "    for dir_type in ['lat', 'lon'] :\n", 
     
    740850  { 
    741851   "cell_type": "code", 
    742    "execution_count": null, 
    743    "id": "6bfcc709", 
    744    "metadata": {}, 
    745    "outputs": [], 
    746    "source": [ 
    747     "# Why is this necessary ?\n", 
    748     "ds['nav_lon_grid_T'].attrs['bounds'] = nav_lon_grid_T.attrs['bounds']\n", 
    749     "ds['nav_lat_grid_T'].attrs['bounds'] = nav_lat_grid_T.attrs['bounds']\n", 
    750     "ds['nav_lon_grid_U'].attrs['bounds'] = nav_lon_grid_U.attrs['bounds']\n", 
    751     "ds['nav_lat_grid_U'].attrs['bounds'] = nav_lat_grid_U.attrs['bounds']\n", 
    752     "ds['nav_lon_grid_V'].attrs['bounds'] = nav_lon_grid_V.attrs['bounds']\n", 
    753     "ds['nav_lat_grid_V'].attrs['bounds'] = nav_lat_grid_V.attrs['bounds']\n", 
    754     "ds['nav_lon_grid_W'].attrs['bounds'] = nav_lon_grid_W.attrs['bounds']\n", 
    755     "ds['nav_lat_grid_W'].attrs['bounds'] = nav_lat_grid_W.attrs['bounds']\n", 
    756     "ds['nav_lon_grid_F'].attrs['bounds'] = nav_lon_grid_F.attrs['bounds']\n", 
    757     "ds['nav_lat_grid_F'].attrs['bounds'] = nav_lat_grid_F.attrs['bounds']" 
    758    ] 
    759   }, 
    760   { 
    761    "cell_type": "code", 
    762    "execution_count": 29, 
     852   "execution_count": 24, 
    763853   "id": "41391c22", 
    764854   "metadata": {}, 
Note: See TracChangeset for help on using the changeset viewer.