Changeset 5939 for TOOLS/MOSAIX
- Timestamp:
- 10/11/21 14:30:35 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
TOOLS/MOSAIX/Build_coordinates_mask.ipynb
r5930 r5939 33 33 { 34 34 "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", 35 60 "id": "26422e3f", 36 61 "metadata": {}, … … 41 66 { 42 67 "cell_type": "code", 43 "execution_count": 2,68 "execution_count": 3, 44 69 "id": "deb66d3d", 45 70 "metadata": {}, … … 47 72 "source": [ 48 73 "#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, 56 83 "id": "2475c8b7", 57 84 "metadata": {}, … … 59 86 "source": [ 60 87 "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", 63 90 " #n_coord = 'https://vesg.ipsl.upmc.fr/thredds/dodsC/work/p86mart/GRAF/DATA/eORCA_R1_coordinates_v1.0.nc'\n", 64 91 " #n_bathy = 'https://vesg.ipsl.upmc.fr/thredds/dodsC/work/p86mart/GRAF/DATA/eORCA_R1_bathy_meter_v2.2.nc'\n", 65 92 " nperio = 6\n", 66 " n_out = ' eORCA_R1_coordinates_mask_test.nc'\n",93 " n_out = 'My_eORCA_R1_coordinates_mask_test.nc'\n", 67 94 "\n", 68 95 "if model == 'ORCA2.3' :\n", 69 96 " 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", 72 99 " n_out = 'My_test.nc'\n", 73 100 " \n", … … 76 103 " n_coord = 'ORCA2.4_coordinates.nc'\n", 77 104 " 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" 79 113 ] 80 114 }, … … 85 119 "source": [ 86 120 "## Open input files\n", 87 "Do not decode time which might be buggy in so emfiles\n",121 "Do not decode time which might be buggy in some files\n", 88 122 "Suppress singleton dimensions if necessary" 89 123 ] … … 91 125 { 92 126 "cell_type": "code", 93 "execution_count": 6,127 "execution_count": 5, 94 128 "id": "0b7c7279", 95 129 "metadata": { … … 104 138 { 105 139 "cell_type": "code", 106 "execution_count": 7,140 "execution_count": 6, 107 141 "id": "9fb92919", 108 142 "metadata": {}, … … 112 146 "output_type": "stream", 113 147 "text": [ 114 " success\n"148 "failed to suppress time\n" 115 149 ] 116 150 } … … 144 178 { 145 179 "cell_type": "code", 146 "execution_count": 8,180 "execution_count": 7, 147 181 "id": "9dfd90e6", 148 182 "metadata": { … … 154 188 "output_type": "stream", 155 189 "text": [ 156 " Exception\n"190 "Normal case\n" 157 191 ] 158 192 } 159 193 ], 160 194 "source": [ 161 "Bathymetry = f_bathy['Bathymetry'] \n",195 "Bathymetry = f_bathy['Bathymetry'].copy()\n", 162 196 "\n", 163 197 "try :\n", … … 176 210 { 177 211 "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", 178 236 "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,201 237 "id": "7e7a49d5", 202 238 "metadata": {}, … … 206 242 "output_type": "stream", 207 243 "text": [ 208 "nperio specified : 4\n"244 "nperio specified : 6\n" 209 245 ] 210 246 } … … 229 265 { 230 266 "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", 231 287 "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" 245 305 ] 246 306 }, … … 248 308 "cell_type": "code", 249 309 "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,268 310 "id": "f0b656a1", 269 311 "metadata": {}, … … 275 317 "mask_W = mask_T\n", 276 318 "\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", 281 323 "\n", 282 324 "mask_U = mask_U.rename ( {'y_grid_T' : 'y_grid_U' , 'x_grid_T' : 'x_grid_U', \n", … … 306 348 { 307 349 "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, 329 351 "id": "274e4321", 330 352 "metadata": {}, 331 353 "outputs": [], 332 354 "source": [ 333 "# Plus compact que ci-dessus, mais un peu incompréhensible .... !!!\n",334 355 "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 335 356 " MaskName = 'mask_' + cd_type\n", 336 357 " 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" 339 361 ] 340 362 }, … … 349 371 { 350 372 "cell_type": "code", 351 "execution_count": 1 6,373 "execution_count": 14, 352 374 "id": "f46e1c93", 353 375 "metadata": {}, … … 377 399 { 378 400 "cell_type": "code", 379 "execution_count": 1 7,401 "execution_count": 15, 380 402 "id": "a7ce1490", 381 403 "metadata": {}, 382 404 "outputs": [], 383 405 "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", 396 428 "\n", 397 429 "nav_lon_grid_T = nav_lon_grid_T.rename( {'y':'y_grid_T', 'x':'x_grid_T'} )\n", … … 406 438 "nav_lat_grid_F = nav_lat_grid_F.rename( {'y':'y_grid_F', 'x':'x_grid_F'} )\n", 407 439 "\n", 440 "# remove _FillValue and missing_value\n", 408 441 "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 409 442 " for dir_type in ['lat', 'lon'] :\n", … … 415 448 { 416 449 "cell_type": "code", 417 "execution_count": 1 8,450 "execution_count": 16, 418 451 "id": "516b4c5c", 419 452 "metadata": {}, … … 442 475 { 443 476 "cell_type": "code", 444 "execution_count": 27,477 "execution_count": 17, 445 478 "id": "de0e6514", 446 479 "metadata": {}, … … 496 529 { 497 530 "cell_type": "code", 498 "execution_count": 20,531 "execution_count": 18, 499 532 "id": "8e4d997c", 500 533 "metadata": {}, 501 534 "outputs": [], 502 535 "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", 503 541 "# remove zero values from areas\n", 504 542 "# need to be define for the extended grid south of -80S\n", 505 543 "# 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, 515 550 "id": "175bae75", 516 551 "metadata": {}, … … 518 553 "source": [ 519 554 "# 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, 535 605 "id": "4f6cfbca", 536 606 "metadata": {}, 537 607 "outputs": [], 538 608 "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", 545 622 "\n", 546 623 "area_grid_T.name = 'area_grid_T'\n", … … 558 635 "area_grid_T.attrs['standard_name'] = 'cell_area'\n", 559 636 "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", 561 638 "area_grid_U.attrs['standard_name'] = 'cell_area'\n", 562 639 "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", 564 641 "area_grid_V.attrs['standard_name'] = 'cell_area'\n", 565 642 "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", 567 644 "area_grid_W.attrs['standard_name'] = 'cell_area'\n", 568 645 "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", 570 647 "area_grid_F.attrs['standard_name'] = 'cell_area'\n", 571 648 "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" 573 654 ] 574 655 }, … … 583 664 { 584 665 "cell_type": "code", 585 "execution_count": 2 3,666 "execution_count": 21, 586 667 "id": "69c404d9", 587 668 "metadata": {}, … … 597 678 " if cdgrd in ['T', 'W'] : \n", 598 679 " 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", 600 682 " if cdgrd == 'U' : \n", 601 683 " 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", 603 686 " if cdgrd == 'V' : \n", 604 687 " 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", 606 690 " if cdgrd == 'F' : \n", 607 691 " 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", 611 697 " dims = ['y_grid_' + cdgrd, 'x_grid_' + cdgrd, 'nvertex_grid_' + cdgrd]\n", 612 698 " \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", 615 701 " \n", 616 702 " idx = [(jcnr,icnr), (jcnr,icnr+1), (jcnr+1,icnr+1), (jcnr+1,icnr)]\n", … … 619 705 " # and complete with periodicity\n", 620 706 " 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", 625 730 " \n", 626 731 " # Rotate cells at the north fold\n", 627 732 " if nperio >= 3 :\n", 628 733 " # Working array for location of northfold\n", 629 " z_fld = nemo.lbc (np.ones ((jpj, jpi)), nperio=nperio, cd_type=cdgrd, psgn=-1. 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", 636 741 " \n", 637 742 " # Invert cells at the symmetric equator\n", 638 743 " 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", 646 751 " bounds_lon.attrs['units'] = 'degrees_east'\n", 647 752 " bounds_lat.attrs['units'] = 'degrees_north'\n", 648 753 " bounds_lon.name = 'bounds_lon_grid_' + cdgrd\n", 649 754 " 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", 650 758 "\n", 651 759 " return bounds_lon, bounds_lat" … … 654 762 { 655 763 "cell_type": "code", 656 "execution_count": 2 4,764 "execution_count": 22, 657 765 "id": "87a8f592", 658 766 "metadata": {}, … … 676 784 { 677 785 "cell_type": "code", 678 "execution_count": 2 8,786 "execution_count": 23, 679 787 "id": "ca87fa1b", 680 788 "metadata": {}, … … 708 816 " 'bounds_lat_grid_F': bounds_lat_grid_F,\n", 709 817 "})\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", 710 820 "for cd_type in ['T', 'U', 'V', 'F', 'W'] :\n", 711 821 " for dir_type in ['lat', 'lon'] :\n", … … 740 850 { 741 851 "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, 763 853 "id": "41391c22", 764 854 "metadata": {},
Note: See TracChangeset
for help on using the changeset viewer.