Changeset 6764
- Timestamp:
- 03/19/24 13:00:07 (9 months ago)
- Location:
- TOOLS/WATER_BUDGET
- Files:
-
- 1 added
- 5 edited
-
ATM_waterbudget.py (modified) (10 diffs)
-
OCE_waterbudget.py (modified) (5 diffs)
-
WaterUtils.py (modified) (9 diffs)
-
libIGCM_date.py (added)
-
libIGCM_sys.py (modified) (3 diffs)
-
lmdz.py (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
TOOLS/WATER_BUDGET/ATM_waterbudget.py
r6688 r6764 50 50 FullIniFile = 'ATM_' + IniFile 51 51 52 print (" Output file : ", FullIniFile )52 print ("FullIniFile : ", FullIniFile ) 53 53 54 54 ## Experiment parameters : read in .ini file 55 55 ## ----------------------------------------- 56 dpar = wu.ReadConfig ( IniFile)56 dpar = wu.ReadConfig (IniFile) 57 57 58 58 ## Configure all needed parameter from existant parameters 59 59 ## ------------------------------------------------------- 60 dpar = wu.SetDatesAndFiles ( dpar)60 dpar = wu.SetDatesAndFiles (dpar) 61 61 62 62 ## Output file with water budget diagnostics … … 77 77 ## Useful functions 78 78 ## ---------------- 79 if repr (readPrec) == "<class 'numpy.float64'>" or readPrec == float :79 if repr (readPrec) == "<class 'numpy.float64'>" or readPrec == float : 80 80 def rprec (ptab) : 81 81 '''This version does nothing … … 121 121 ## Open history files 122 122 ## ------------------ 123 d_ATM_his = xr.open_dataset ( file_ATM_his, use_cftime=True, decode_times=True, decode_cf=True ) .squeeze()124 if SECHIBA : 125 d_SRF_his = xr.open_dataset ( file_SRF_his, use_cftime=True, decode_times=True, decode_cf=True ) .squeeze()123 d_ATM_his = xr.open_dataset ( file_ATM_his, use_cftime=True, decode_times=True, decode_cf=True )#.squeeze() 124 if SECHIBA : 125 d_SRF_his = xr.open_dataset ( file_SRF_his, use_cftime=True, decode_times=True, decode_cf=True )#.squeeze() 126 126 if Routing == 'SECHIBA' : d_RUN_his = d_SRF_his 127 if Routing == 'SIMPLE' : d_RUN_his = xr.open_dataset ( file_RUN_his, use_cftime=True, decode_times=True, decode_cf=True ) .squeeze()127 if Routing == 'SIMPLE' : d_RUN_his = xr.open_dataset ( file_RUN_his, use_cftime=True, decode_times=True, decode_cf=True )#.squeeze() 128 128 129 129 echo ( f'{file_ATM_his = }' ) … … 136 136 echo ( f'\nRun length : {(dtime/np.timedelta64(1, "D")).values:8.2f} days' ) 137 137 dtime_sec = (dtime/np.timedelta64(1, "s")).values.item() # Convert in seconds 138 dpar['Experiment']['dtime_sec'] = dtime_sec139 138 140 139 ##-- Compute length of each period 141 dtime_per = (d_ATM_his.time_counter_bounds[ :,-1] - d_ATM_his.time_counter_bounds[:,0] )140 dtime_per = (d_ATM_his.time_counter_bounds[...,-1] - d_ATM_his.time_counter_bounds[...,0] ) 142 141 echo ( '\nPeriods lengths (days) : {:} days'.format ( (dtime_per/np.timedelta64(1, "D")).values ) ) 143 142 dtime_per_sec = (dtime_per/np.timedelta64(1, "s")).values # In seconds 144 dtime_per_sec = xr.DataArray (dtime_per_sec, dims= ["time_counter", ], coords=[d_ATM_his.time_counter,])143 dtime_per_sec = xr.DataArray (dtime_per_sec, dims=("time_counter",), coords=(d_ATM_his.time_counter,) ) 145 144 dtime_per_sec.attrs['unit'] = 's' 146 dpar['Experiment']['dtime_per_sec'] = dtime_per_sec147 145 148 146 ##-- Number of years (approximative) … … 393 391 394 392 echo ( '\n====================================================================================' ) 393 echo ( f'-- ATM checking GRAV -- {Title} ' ) 394 ## Estimation de GRAV 395 zGrav = d_ATM_his.psol / d_ATM_his.mass.sum(dim='presnivs') 396 echo ( f'Grav estimated {zGrav.min().values} {zGrav.max().values} versus specified {GRAV}' ) 397 398 if GRAV > zGrav.max() or GRAV < zGrav.min() : 399 raise ValueError ( 'Error of GRAV value' ) 400 401 echo ( '\n====================================================================================' ) 395 402 echo ( f'-- ATM changes in stores -- {Title} ' ) 396 403 … … 401 408 ATM.Ahyb = d_ATM_his['Ahyb'].squeeze() 402 409 ATM.Bhyb = d_ATM_his['Bhyb'].squeeze() 410 411 #DYN.Ahyb_beg = d_DYN_beg['ap'] 412 #DYN.Bhyb_beg = d_DYN_beg['bp'] 413 #DYN.Ahyb_end = d_DYN_end['ap'] 414 #DYN.Bhyb_end = d_DYN_end['bp'] 403 415 404 416 echo ( 'Surface pressure' ) … … 413 425 DYN.pres_beg = ATM.Ahyb + ATM.Bhyb * DYN.psol_beg 414 426 DYN.pres_end = ATM.Ahyb + ATM.Bhyb * DYN.psol_end 427 428 #DYN.pres_beg = DYN.Ahyb_beg + DYN.Bhyb_beg * DYN.psol_beg 429 #DYN.pres_end = DYN.Ahyb_end + DYN.Bhyb_end * DYN.psol_end 415 430 416 431 echo ( 'Check computation of pressure levels' ) … … 427 442 ind[7] = (DYN.pres_end[-1]).max().item() 428 443 429 if any ( 444 if any (ind != 0) : 430 445 echo ( 'All values should be zero' ) 431 446 echo ( f'(DYN.pres_beg[ 0]-DYN.psol_beg).min().item() = {ind[0]}' ) … … 457 472 458 473 echo ( 'Atmosphere mass (kg)' ) 459 DYN.mass_air_beg = DYN_stock_int ( DYN.mass_beg_2D ) 460 DYN.mass_air_end = DYN_stock_int ( DYN.mass_end_2D ) 474 DYN.mass_air_beg = DYN_stock_int ( DYN.mass_beg_2D ) 475 DYN.mass_air_end = DYN_stock_int ( DYN.mass_end_2D ) 476 477 DYN.mass_psol_beg = DYN_stock_int ( DYN.psol_beg ) / GRAV 478 DYN.mass_psol_end = DYN_stock_int ( DYN.psol_end ) / GRAV 461 479 462 480 echo ( 'Atmosphere mass change (kg)' ) 463 DYN.diff_mass_air = DYN.mass_air_end - DYN.mass_air_beg 481 DYN.diff_mass_air = DYN.mass_air_end - DYN.mass_air_beg 482 DYN.diff_mass_psol = DYN.mass_psol_end - DYN.mass_psol_beg 483 484 echo ( f'\nChk psol vs. mass (dynamics) -- {Title} ' ) 485 echo ( '------------------------------------------------------------------------------------' ) 486 echo ( 'DYN.mass_air_beg = {:12.6e} kg | DYN.mass_psol_beg = {:12.6e} kg | diff: {:12.6e}'.format (DYN.mass_air_beg, DYN.mass_psol_beg, DYN.mass_air_beg-DYN.mass_psol_beg) ) 487 echo ( 'DYN.mass_air_end = {:12.6e} kg | DYN.mass_psol_end = {:12.6e} kg | diff: {:12.6e}'.format (DYN.mass_air_end, DYN.mass_psol_end, DYN.mass_air_end-DYN.mass_psol_end) ) 488 echo ( 'DYN.diff_mass_air = {:12.6e} kg'.format (DYN.diff_mass_air ) ) 489 echo ( 'DYN.diff_mass_psol = {:12.6e} kg'.format (DYN.diff_mass_psol) ) 490 echo ( 'diff = {:12.6e} kg'.format (DYN.diff_mass_psol-DYN.diff_mass_air) ) 464 491 465 492 echo ( f'\nChange of atmosphere mass (dynamics) -- {Title} ' ) … … 469 496 echo ( f'dMass (atm) : {(DYN.diff_mass_air/DYN.mass_air_beg):9.4e} (-)' ) 470 497 498 echo ( f'\nChange of atmosphere mass from psol (dynamics) -- {Title} ' ) 499 echo ( '------------------------------------------------------------------------------------' ) 500 echo ( 'DYN.mass_psol_beg = {:12.6e} kg | DYN.mass_psol_end = {:12.6e} kg'.format (DYN.mass_psol_beg, DYN.mass_psol_end) ) 501 echo ( f'dMass (atm) : {DYN.diff_mass_psol:9.4e} kg' ) 502 echo ( f'dMass (atm) : {(DYN.diff_mass_psol/DYN.mass_psol_beg):9.4e} (-)' ) 503 471 504 echo ( 'Vertical and horizontal integral, and sum of liquid, solid and vapor water phases' ) 472 505 if LMDZ : 473 506 if 'H2Ov' in d_DYN_beg.variables : 474 507 echo ('reading LATLON : H2Ov, H2Ol, H2Oi' ) 475 DYN.wat_beg = lmdz.geo3point ( (d_DYN_beg['H2Ov'] + d_DYN_beg['H2Ol'] + d_DYN_beg['H2Oi']).isel(rlonv=slice(0,-1) ), dim1d='cell' ) 476 DYN.wat_end = lmdz.geo3point ( (d_DYN_end['H2Ov'] + d_DYN_end['H2Ol'] + d_DYN_end['H2Oi']).isel(rlonv=slice(0,-1) ), dim1d='cell' ) 508 DYN.wat_beg_gas = lmdz.geo2point ( d_DYN_beg['H2Ov'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 509 DYN.wat_beg_liq = lmdz.geo2point ( d_DYN_beg['H2Ol'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 510 DYN.wat_beg_sol = lmdz.geo2point ( d_DYN_beg['H2Oi'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 511 DYN.wat_end_gas = lmdz.geo2point ( d_DYN_end['H2Ov'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 512 DYN.wat_end_liq = lmdz.geo2point ( d_DYN_end['H2Ol'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 513 DYN.wat_end_sol = lmdz.geo2point ( d_DYN_end['H2Oi'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 477 514 if 'H2O_g' in d_DYN_beg.variables : 478 515 echo ('reading LATLON : H2O_g, H2O_l, H2O_s' ) 479 DYN.wat_beg = lmdz.geo3point ( (d_DYN_beg['H2O_g'] + d_DYN_beg['H2O_l'] + d_DYN_beg['H2O_s']).isel(rlonv=slice(0,-1) ), dim1d='cell' ) 480 DYN.wat_end = lmdz.geo3point ( (d_DYN_end['H2O_g'] + d_DYN_end['H2O_l'] + d_DYN_end['H2O_s']).isel(rlonv=slice(0,-1) ), dim1d='cell' ) 516 DYN.wat_beg_gas = lmdz.geo2point ( d_DYN_beg['H2O_g'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 517 DYN.wat_beg_liq = lmdz.geo2point ( d_DYN_beg['H2O_l'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 518 DYN.wat_beg_sol = lmdz.geo2point ( d_DYN_beg['H2O_s'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 519 DYN.wat_end_gas = lmdz.geo2point ( d_DYN_end['H2O_g'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 520 DYN.wat_end_liq = lmdz.geo2point ( d_DYN_end['H2O_l'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 521 DYN.wat_end_sol = lmdz.geo2point ( d_DYN_end['H2O_i'].isel(rlonv=slice(0,-1) ), dim1d='cell' ) 481 522 if ICO : 482 if 523 if 'H2Ov_g' in d_DYN_beg.variables : 483 524 echo ('reading ICO : H2Ov_g, H2Ov_l, H2Ov_s' ) 484 DYN.wat_beg = d_DYN_beg['H2Ov_g'] + d_DYN_beg['H2Ov_l'] + d_DYN_beg['H2Ov_s'] 485 DYN.wat_end = d_DYN_end['H2Ov_g'] + d_DYN_end['H2Ov_l'] + d_DYN_end['H2Ov_s'] 525 DYN.wat_beg_gas = d_DYN_beg['H2Ov_g'] 526 DYN.wat_beg_liq = d_DYN_beg['H2Ov_l'] 527 DYN.wat_beg_sol = d_DYN_beg['H2Ov_s'] 528 DYN.wat_end_gas = d_DYN_end['H2Ov_g'] 529 DYN.wat_end_liq = d_DYN_end['H2Ov_l'] 530 DYN.wat_end_sol = d_DYN_end['H2Ov_s'] 486 531 if 'H2O_g' in d_DYN_beg.variables : 487 532 echo ('reading ICO : H2O_g, H2O_l, H2O_s' ) 488 DYN.wat_beg = d_DYN_beg['H2O_g'] + d_DYN_beg['H2O_l'] + d_DYN_beg['H2O_s'] 489 DYN.wat_end = d_DYN_end['H2O_g'] + d_DYN_end['H2O_l'] + d_DYN_end['H2O_s'] 490 if 'q' in d_DYN_beg.variables : 491 echo ('reading ICO : q' ) 492 DYN.wat_beg = d_DYN_beg['q'].isel(nq=0) + d_DYN_beg['q'].isel(nq=1) + d_DYN_beg['q'].isel(nq=2) 493 DYN.wat_end = d_DYN_end['q'].isel(nq=0) + d_DYN_end['q'].isel(nq=1) + d_DYN_end['q'].isel(nq=2) 494 495 if 'lev' in DYN.wat_beg.dims : 496 DYN.wat_beg = DYN.wat_beg.rename ( {'lev':'sigs'} ) 497 DYN.wat_end = DYN.wat_end.rename ( {'lev':'sigs'} ) 498 533 DYN.wat_beg_gas = d_DYN_beg['H2O_g'] 534 DYN.wat_beg_liq = d_DYN_beg['H2O_l'] 535 DYN.wat_beg_sol = d_DYN_beg['H2O_s'] 536 DYN.wat_end_gas = d_DYN_end['H2O_g'] 537 DYN.wat_end_liq = d_DYN_end['H2O_l'] 538 DYN.wat_end_sol = d_DYN_end['H2O_s'] 539 #if 'q' in d_DYN_beg.variables : 540 # echo ('reading ICO : q' ) 541 # DYN.wat_beg_gas = d_DYN_beg['q'].isel(nq=0) 542 # DYN.wat_beg_liq = d_DYN_beg['q'].isel(nq=1) 543 # DYN.wat_beg_sol = d_DYN_beg['q'].isel(nq=2) 544 # DYN.wat_end_gas = d_DYN_end['q'].isel(nq=0) 545 # DYN.wat_end_liq = d_DYN_end['q'].isel(nq=1) 546 # DYN.wat_end_sol = d_DYN_end['q'].isel(nq=2) 547 548 if 'lev' in DYN.wat_beg_gas.dims : 549 DYN.wat_beg_gas = DYN.wat_beg_gas.rename ( {'lev':'sigs'} ) 550 DYN.wat_beg_liq = DYN.wat_beg_liq.rename ( {'lev':'sigs'} ) 551 DYN.wat_beg_sol = DYN.wat_beg_sol.rename ( {'lev':'sigs'} ) 552 DYN.wat_end_gas = DYN.wat_end_gas.rename ( {'lev':'sigs'} ) 553 DYN.wat_end_liq = DYN.wat_end_liq.rename ( {'lev':'sigs'} ) 554 DYN.wat_end_sol = DYN.wat_end_sol.rename ( {'lev':'sigs'} ) 555 556 DYN.wat_beg = DYN.wat_beg_gas + DYN.wat_beg_liq + DYN.wat_beg_sol 557 DYN.wat_end = DYN.wat_end_gas + DYN.wat_end_liq + DYN.wat_end_sol 558 499 559 echo ( 'Compute water content : vertical and horizontal integral' ) 500 501 DYN.wat_beg_2D = (DYN.mass_beg * DYN.wat_beg).sum (dim='sigs') 502 DYN.wat_end_2D = (DYN.mass_end * DYN.wat_end).sum (dim='sigs') 503 504 DYN.mass_wat_beg = DYN_stock_int ( DYN.wat_beg_2D ) 505 DYN.mass_wat_end = DYN_stock_int ( DYN.wat_end_2D ) 560 DYN.wat_beg_gas_2D = (DYN.mass_beg * DYN.wat_beg_gas).sum (dim='sigs') 561 DYN.wat_beg_liq_2D = (DYN.mass_beg * DYN.wat_beg_liq).sum (dim='sigs') 562 DYN.wat_beg_sol_2D = (DYN.mass_beg * DYN.wat_beg_sol).sum (dim='sigs') 563 DYN.wat_end_gas_2D = (DYN.mass_end * DYN.wat_end_gas).sum (dim='sigs') 564 DYN.wat_end_liq_2D = (DYN.mass_end * DYN.wat_end_liq).sum (dim='sigs') 565 DYN.wat_end_sol_2D = (DYN.mass_end * DYN.wat_end_sol).sum (dim='sigs') 566 DYN.wat_beg_2D = (DYN.mass_beg * DYN.wat_beg ).sum (dim='sigs') 567 DYN.wat_end_2D = (DYN.mass_end * DYN.wat_end ).sum (dim='sigs') 568 569 DYN.mass_wat_beg_gas = DYN_stock_int ( DYN.wat_beg_gas_2D ) 570 DYN.mass_wat_beg_liq = DYN_stock_int ( DYN.wat_beg_liq_2D ) 571 DYN.mass_wat_beg_sol = DYN_stock_int ( DYN.wat_beg_sol_2D ) 572 DYN.mass_wat_end_gas = DYN_stock_int ( DYN.wat_end_gas_2D ) 573 DYN.mass_wat_end_liq = DYN_stock_int ( DYN.wat_end_liq_2D ) 574 DYN.mass_wat_end_sol = DYN_stock_int ( DYN.wat_end_sol_2D ) 575 DYN.mass_wat_beg = DYN_stock_int ( DYN.wat_beg_2D ) 576 DYN.mass_wat_end = DYN_stock_int ( DYN.wat_end_2D ) 506 577 507 578 echo ( 'Variation of water content' ) 508 DYN.diff_mass_wat = DYN.mass_wat_end - DYN.mass_wat_beg 509 579 DYN.diff_mass_wat_gas = DYN.mass_wat_end_gas - DYN.mass_wat_beg_gas 580 DYN.diff_mass_wat_liq = DYN.mass_wat_end_liq - DYN.mass_wat_beg_liq 581 DYN.diff_mass_wat_sol = DYN.mass_wat_end_sol - DYN.mass_wat_beg_sol 582 DYN.diff_mass_wat = DYN.mass_wat_end - DYN.mass_wat_beg 583 584 if 'prw_end' in d_ATM_his : 585 echo ( 'Check with pr variables' ) 586 if ATM_HIS == 'latlon' : 587 echo ( ' latlon case' ) 588 ATM.wat_prw = lmdz.geo2point ( rprec (d_ATM_his ['prw_end'] [-1]), dim1d='cell' ) 589 ATM.wat_prlw = lmdz.geo2point ( rprec (d_ATM_his ['prlw_end'] [-1]), dim1d='cell' ) 590 ATM.wat_prsw = lmdz.geo2point ( rprec (d_ATM_his ['prsw_end'] [-1]), dim1d='cell' ) 591 ATM.wat_prbsw = lmdz.geo2point ( rprec (d_ATM_his ['prbsw_end'][-1]), dim1d='cell' ) 592 if ATM_HIS == 'ico' : 593 echo ( ' ico case' ) 594 ATM.wat_prw = rprec (d_ATM_his ['prw_end'] [-1]) 595 ATM.wat_prlw = rprec (d_ATM_his ['prlw_end'] [-1]) 596 ATM.wat_prsw = rprec (d_ATM_his ['prsw_end'] [-1]) 597 ATM.wat_prbsw = rprec (d_ATM_his ['prbsw_end'][-1]) 598 599 ATM.mass_wat_prw = DYN_stock_int ( ATM.wat_prw ) 600 ATM.mass_wat_prlw = DYN_stock_int ( ATM.wat_prlw ) 601 ATM.mass_wat_prsw = DYN_stock_int ( ATM.wat_prsw ) 602 ATM.mass_wat_prbsw = DYN_stock_int ( ATM.wat_prbsw ) 603 604 ATM.mass_wat_pr = ATM.mass_wat_prw + ATM.mass_wat_prlw + ATM.mass_wat_prsw + ATM.mass_wat_prbsw 605 606 else : 607 ATM.wat_pr = 0 608 ATM.wat_prw = 0 609 ATM.wat_prlw = 0 610 ATM.wat_prsw = 0 611 ATM.wat_prbsw = 0 612 ATM.mass_wat_prw = 0 613 ATM.mass_wat_prlw = 0 614 ATM.mass_wat_prsw = 0 615 ATM.mass_wat_prbsw = 0 616 ATM.mass_wat_pr = 0 617 618 echo ( f'\nCheck atmosphere water content -- {Title} ' ) 619 echo ( '------------------------------------------------------------------------------------' ) 620 echo ( 'DYN.mass_wat_end_gas = {:12.6e} kg | ATM.mass_wat_prw = {:12.6e} kg'.format (DYN.mass_wat_end_gas, ATM.mass_wat_prw ) ) 621 echo ( 'DYN.mass_wat_end_liq = {:12.6e} kg | ATM.mass_wat_prlw = {:12.6e} kg'.format (DYN.mass_wat_end_liq, ATM.mass_wat_prlw ) ) 622 echo ( 'DYN.mass_wat_end_sol = {:12.6e} kg | ATM.mass_wat_prsw = {:12.6e} kg'.format (DYN.mass_wat_end_sol, ATM.mass_wat_prsw ) ) 623 echo ( 'DYN.mass_wat_end_sno = {:12.6e} kg | ATM.mass_wat_prbsw = {:12.6e} kg'.format ( 0. , ATM.mass_wat_prbsw) ) 624 echo ( 'DYN.mass_wat_end = {:12.6e} kg | ATM.mass_wat_pr = {:12.6e} kg'.format (DYN.mass_wat_end , ATM.mass_wat_pr ) ) 625 prtFlux ( 'Error water mass (dyn - atm)', DYN.mass_wat_end - ATM.mass_wat_pr, 'e', True ) 626 510 627 echo ( f'\nChange of atmosphere water content (dynamics) -- {Title} ' ) 511 628 echo ( '------------------------------------------------------------------------------------' ) -
TOOLS/WATER_BUDGET/OCE_waterbudget.py
r6688 r6764 1 #!/usr/bin/env python31 !/usr/bin/env python3 2 2 ### 3 3 ### Script to check water conservation in NEMO … … 51 51 52 52 53 print ("Output file : ", FullIniFile)53 print ("Output parameter file : ", FullIniFile) 54 54 55 55 ## Experiment parameters … … 59 59 ## Configure all needed parameter from existant parameters 60 60 ## ------------------------------------------------------- 61 dpar = wu.SetDatesAndFiles ( 61 dpar = wu.SetDatesAndFiles (dpar) 62 62 63 63 ## Output file with water budget diagnostics … … 80 80 81 81 # Degrades precision 82 if repr (readPrec) == "<class 'numpy.float64'>" or readPrec == float :82 if repr (readPrec) == "<class 'numpy.float64'>" or readPrec == float : 83 83 def rprec (ptab) : 84 84 '''This version does nothing … … 271 271 272 272 ErrorCount = 0 273 274 273 ErrorCount += extract_and_rebuild ( file_name=file_OCE_beg, tar_restart=tar_restart_beg, file_dir_comp=FileDirOCE ) 275 274 ErrorCount += extract_and_rebuild ( file_name=file_OCE_end, tar_restart=tar_restart_end, file_dir_comp=FileDirOCE ) -
TOOLS/WATER_BUDGET/WaterUtils.py
r6688 r6764 31 31 import lmdz 32 32 33 def ReadConfig ( ini_file, default_ini_file='defaults.ini') :33 def ReadConfig (ini_file, default_ini_file='defaults.ini') : 34 34 ''' 35 35 Reads experiment parameters … … 92 92 return config2dict (params) 93 93 94 def SetDatesAndFiles ( pdict) :94 def SetDatesAndFiles (pdict) : 95 95 ''' 96 96 From readed experiment parameters, set all needed parameters … … 126 126 127 127 pdict['Files']['TmpDir'] = mm['TmpDir'] 128 pdict['libIGCM'].update (mm )128 pdict['libIGCM'].update (mm.__dict__) 129 129 130 130 ## Complete configuration … … 156 156 pdict['Config']['readPrec'] = np.float16 157 157 158 159 158 ## Defines begining and end of experiment 160 159 ## -------------------------------------- … … 162 161 pdict['Experiment']['DateBegin'] = f"{pdict['Experiment']['YearBegin']}0101" 163 162 else : 164 YearBegin, MonthBegin, DayBegin = ldate.GetYearMonthDay ( pdict['Experiment']['DateBegin'] ) 165 DateBegin = FormatToGregorian (pdict['Experiment']['DateBegin']) 166 pdict['Experiment']['YearBegin'] = YearBegin 163 pdict['Experiment']['YearBegin'],pdict['Experiment']['MonthBegin'], pdict['Experiment']['DayBegin'] = ldate.GetYearMonthDay ( pdict['Experiment']['DateBegin'] ) 164 pdict['Experiment']['DateBegin'] = ldate.ConvertFormatToGregorian (pdict['Experiment']['DateBegin']) 167 165 168 166 if not pdict['Experiment']['DateEnd'] : 169 167 pdict['Experiment']['DateEnd'] = f"{pdict['Experiment']['YearEnd']}1231" 170 168 else : 171 YearEnd, MonthEnd, DayEnd = ldate.GetYearMonthDay ( pdict['Experiment']['DateEnd'] ) 172 DateEnd = FormatToGregorian (pdict['Experiment']['DateEnd']) 173 pdict['Experiment']['DateEnd'] = DateEnd 169 pdict['Experiment']['YearEnd'], pdict['Experiment']['MonthEnd'], pdict['Experiment']['DayEnd'] = ldate.GetYearMonthDay ( pdict['Experiment']['DateEnd'] ) 170 pdict['Experiment']['DateEnd'] = ldate.ConvertFormatToGregorian (pdict['Experiment']['DateEnd']) 174 171 175 172 if not pdict['Experiment']['PackFrequency'] : 176 PackFrequencypdict['Experiment']['PackFrequency'] = f"{YearEnd - YearBegin+ 1}YE"173 pdict['Experiment']['PackFrequency'] = f"{pdict['Experiment']['YearEnd'] - pdict['Experiment']['YearBegin'] + 1}YE" 177 174 178 179 175 ## Set directory to extract files 180 176 ## ------------------------------ 181 177 if not pdict['Files']['FileDir'] : 182 178 pdict['Files']['FileDir'] = os.path.join ( pdict['Files']['TmpDir'], f"WATER_{pdict['Experiment']['JobName']}" ) 183 184 179 if not os.path.isdir ( pdict['Files']['FileDir'] ) : os.makedirs ( pdict['Files']['FileDir'] ) 185 180 181 if not pdict['Files']['FileDirOCE'] : 182 pdict['Files']['FileDirOCE'] = os.path.join ( pdict['Files']['FileDir'], "OCE" ) 183 if not os.path.isdir ( pdict['Files']['FileDirOCE'] ) : os.makedirs ( pdict['Files']['FileDirOCE'] ) 184 185 if not pdict['Files']['FileDirICE'] : 186 pdict['Files']['FileDirICE'] = os.path.join ( pdict['Files']['FileDir'], "ICE" ) 187 if not os.path.isdir ( pdict['Files']['FileDirICE'] ) : os.makedirs ( pdict['Files']['FileDirICE'] ) 188 186 189 FileOut = str2value (pdict['Files']['FileOut']) 187 190 if not FileOut : 188 FileOut = f"ATM_waterbudget_{pdict['Experiment']['JobName']}_{pdict['Experiment'][' YearBegin']}_{pdict['Experiment']['YearEnd']}"191 FileOut = f"ATM_waterbudget_{pdict['Experiment']['JobName']}_{pdict['Experiment']['DateBegin']}_{pdict['Experiment']['DateEnd']}" 189 192 if pdict['Experiment']['ICO'] : 190 193 if pdict['Experiment']['ATM_HIS'] == 'latlon' : FileOut = f'{FileOut}_LATLON' … … 203 206 echo ( f"FileDir : {pdict['Files']['FileDir']}" , f_out=f_out ) 204 207 205 echo ( f"\nDealing with {pdict['libIGCM']['L_EXP'] }", f_out )208 echo ( f"\nDealing with {pdict['libIGCM']['L_EXP']=}", f_out ) 206 209 207 210 ## Creates model output directory names … … 280 283 return pdict 281 284 282 def SetRestartNames ( 285 def SetRestartNames (pdict, f_out) : 283 286 ''' 284 287 Defines dates for restart files … … 309 312 pdict['Files']['TarRestartPeriod_end_DateBeg'] = ldate.DateAddPeriod ( pdict['Files']['TarRestartPeriod_end_DateEnd'], 310 313 f"-{pdict['Experiment']['PackFrequency']}" ) 311 pdict['Files']['TarRestartPeriod_end_DateBeg'] = ldate.AddOneDayToDate ( pdict['Files']['TarRestartPeriod_end_DateBeg'], 312 pdict['Experiment']['CalendarType'] ) 314 pdict['Files']['TarRestartPeriod_end_DateBeg'] = ldate.AddOneDayToDate ( pdict['Files']['TarRestartPeriod_end_DateBeg'], pdict['Experiment']['CalendarType'] ) 313 315 314 316 pdict['Files']['TarRestartPeriod_end'] = f"{pdict['Files']['TarRestartPeriod_end_DateBeg']}_{pdict['Files']['TarRestartPeriod_end_DateEnd']}" … … 624 626 The dictionary must have two levels (see configparser) 625 627 ''' 626 zconf = configparser.ConfigParser ( interpolation=configparser.ExtendedInterpolation())628 zconf = configparser.ConfigParser (interpolation=configparser.ExtendedInterpolation ()) 627 629 zconf.optionxform = str # To keep capitals 628 630 zconf.read_dict(dict2str(pdict)) -
TOOLS/WATER_BUDGET/libIGCM_sys.py
r6688 r6764 8 8 environment variables and commands. 9 9 All those definitions depend on host particularities. 10 It manages a stack mechanism and test validity of operations.11 10 12 11 This software is governed by the CeCILL license under French law and … … 47 46 # return unDefined 48 47 49 def Mach ( 48 def Mach (long=False) : 50 49 ''' 51 50 Find the computer we are on … … 84 83 return zmach 85 84 86 def config ( JobName=None , TagName=None , SpaceName=None, ExperimentName=None, 87 LongName=None, ModelName=None, ShortName=None, 88 Source=None , Host=None , ConfigCard=None, User=None, Group=None, 89 TGCC_User='p86mart', TGCC_Group='gen12006', IDRIS_User='rces009', IDRIS_Group='ces', 90 ARCHIVE=None, SCRATCHDIR=None, STORAGE=None, R_IN=None, R_OUT=None, 91 R_FIG=None, L_EXP=None, 92 R_SAVE=None , R_FIGR=None, R_BUF=None , R_BUFR=None, R_BUF_KSH=None, 93 REBUILD_DIR=None, POST_DIR=None, 94 ThreddsPrefix=None, R_GRAF=None, DB=None, 95 IGCM_OUT=None, IGCM_OUT_name=None, rebuild=None, TmpDir=None, 96 Out='dict') : 97 ''' 98 Defines the libIGCM directories - Returns a dictionnary or a namespace 99 100 Source : for Spip 101 Possibilities : 102 Local : local (~/Data/IGCMG/...) 103 TGCC_sshfs : TGCC disks mounted via sshfs 104 TGCC_thredds : thredds TGCC via IPSL 105 ('https://thredds-su.ipsl.fr/thredds/dodsC/tgcc_thredds/store/...) 106 107 Exemple of use : 108 libIGCM = libIGCM_sys.config ( .... ) 109 globals().update ( libIGCM ) 110 111 ''' 112 #if not IGCM_OUT_name : 113 # if Source == 'TGCC_thredds' : IGCM_OUT_name = '' 114 # else : IGCM_OUT_name = 'IGCM_OUT' 115 116 if not Host : Host = Mach (long=False) 117 118 LocalUser = os.environ ['USER'] 119 120 # =========================================================================================== 121 # Reads config.card if available 122 if ConfigCard : 123 # Text existence of ConfigCard 124 if not os.path.exists (ConfigCard ) : 125 raise FileExistsError ( f"File not found : {ConfigCard = }" ) 126 127 ## Creates parser for reading .ini input file 128 MyReader = configparser.ConfigParser (interpolation=configparser.ExtendedInterpolation() ) 129 MyReader.optionxform = str # To keep capitals 130 131 MyReader.read (ConfigCard) 132 133 JobName = MyReader['UserChoices']['JobName'] 134 #----- Short Name of Experiment 135 ExperimentName = MyReader['UserChoices']['ExperimentName'] 136 #----- DEVT TEST PROD 137 SpaceName = MyReader['UserChoices']['SpaceName'] 138 LongName = MyReader['UserChoices']['LongName'] 139 ModelName = MyReader['UserChoices']['ModelName'] 140 TagName = MyReader['UserChoices']['TagName'] 141 142 ### =========================================================================================== 143 ## Part specific to access by OpenDAP/Thredds server 144 145 # =========================================================================================== 146 if Source == 'TGCC_thredds' : 147 if not User and TGCC_User : User = TGCC_User 148 if not Group and TGCC_Group : Group = TGCC_Group 149 150 if not ThreddsPrefix : 151 ThreddsPrefix = 'https://thredds-su.ipsl.fr/thredds/dodsC/tgcc_thredds' 152 153 if not ARCHIVE : ARCHIVE = f'{ThreddsPrefix}/store/{TGCC_User}' 154 if not R_FIG : R_FIG = f'{ThreddsPrefix}/work/{TGCC_User}' 155 if not R_IN : R_IN = f'{ThreddsPrefix}/work/igcmg/IGCM' 156 if not R_GRAF : R_GRAF = f'{ThreddsPrefix}/work/p86mart/GRAF/DATA' 157 158 # =========================================================================================== 159 if Source == 'IDRIS_thredds' : 160 if not User and IDRIS_User : User = IDRIS_User 161 if not Group and IDRIS_Group : Group = IDRIS_Group 162 163 if not ThreddsPrefix : 164 ThreddsPrefix = 'https://thredds-su.ipsl.fr/thredds/catalog/idris_thredds' 165 166 if not ARCHIVE : ARCHIVE = f'{ThreddsPrefix}/store/{IDRIS_User}' 167 if not R_FIG : R_FIG = f'{ThreddsPrefix}/work/{IDRIS_User}' 168 if not R_IN : R_IN = f'{ThreddsPrefix}/work/igcmg/IGCM' 169 if not R_GRAF : R_GRAF = \ 170 'https://thredds-su.ipsl.fr/thredds/dodsC/tgcc_thredds/work/p86mart/GRAF/DATA' 171 172 ### =========================================================================================== 173 ## Machine dependant part 174 175 # =========================================================================================== 176 if Host == 'Spip' : 177 if not User : User = LocalUser 178 IGCM_OUT_name = 'IGCM_OUT' 179 if not ARCHIVE : ARCHIVE = os.path.join ( Path.home (), 'Data' ) 180 if not SCRATCHDIR : SCRATCHDIR = os.path.join ( Path.home (), 'Scratch' ) 181 if not R_BUF : R_BUF = os.path.join ( Path.home (), 'Scratch' ) 182 if not R_FIG : R_FIG = os.path.join ( Path.home (), 'Data' ) 183 184 if not STORAGE : STORAGE = ARCHIVE 185 if not R_IN : R_IN = os.path.join ( '/', 'Users', 'marti', 'Data', 'IGCM' ) 186 if not R_GRAF : R_GRAF = os.path.join ( '/', 'Users', 'marti', 'GRAF', 'DATA' ) 187 if not DB : DB = os.path.join ( '/', 'Users', 'marti', 'GRAF', 'DB' ) 188 if not TmpDir : TmpDir = os.path.join ( '/', 'Users', LocalUser, 'Scratch' ) 189 190 # =========================================================================================== 191 if ( 'Irene' in Host ) or ( 'Rome' in Host ) : 192 193 LocalHome = subprocess.getoutput ( 'ccc_home --ccchome' ) 194 LocalGroup = os.path.basename ( os.path.dirname (LocalHome)) 195 if not User or User == 'marti' : 196 if not TGCC_User : User = LocalUser 197 else : User = TGCC_User 198 199 if not Group : 200 if TGCC_Group : Group = TGCC_Group 201 else : Group = LocalGroup 202 203 IGCM_OUT_name = 'IGCM_OUT' 204 205 if not R_IN : 206 R_IN = os.path.join ( subprocess.getoutput ( 207 'ccc_home --cccwork -d igcmg -u igcmg' ), 'IGCM') 208 if not ARCHIVE : 209 ARCHIVE = subprocess.getoutput ( 210 f'ccc_home --cccstore -u {User} -d {Group}' ) 211 if not STORAGE : 212 STORAGE = subprocess.getoutput ( 213 f'ccc_home --cccwork -u {User} -d {Group}' ) 214 if not SCRATCHDIR : 215 SCRATCHDIR = subprocess.getoutput ( 216 f'ccc_home --cccscratch -u {User} -d {Group}' ) 217 if not R_BUF : 218 R_BUF = subprocess.getoutput ( 219 f'ccc_home --cccscratch -u {User} -d {Group}' ) 220 if not R_FIG : 221 R_FIG = subprocess.getoutput ( 222 f'ccc_home --cccwork -u {User} -d {Group}' ) 223 if not R_GRAF : 224 R_GRAF = os.path.join ( subprocess.getoutput ( 225 'ccc_home --cccwork -d drf -u p86mart'), 'GRAF', 'DATA' ) 226 if not DB : 227 DB = os.path.join ( subprocess.getoutput ( 228 'ccc_home --cccwork -d igcmg -u igcmg'), 'database' ) 229 if not rebuild : 85 class config : 86 def __str__ (self) : 87 return str (self.__dict__) 88 89 def __repr__ (self) : 90 return str (self.__dict__) 91 92 def __name__ (self) : 93 return self.__class__.__name__ 94 95 def __getitem__ (self, key): 96 return getattr (self, key) 97 98 def __iter__ (self) : 99 return self 100 101 def __next__ (self) : 102 if self.index == 0: 103 raise StopIteration 104 self.index = self.index - 1 105 return self.data[self.index] 106 107 def __init__ (self, JobName=None, TagName=None, SpaceName=None, ExperimentName=None, 108 LongName=None, ModelName=None, ShortName=None, 109 Source=None, Host=None, ConfigCard=None, User=None, Group=None, 110 TGCC_User='p86mart', TGCC_Group='gen12006', IDRIS_User='rces009', IDRIS_Group='ces', 111 ARCHIVE=None, SCRATCHDIR=None, STORAGE=None, R_IN=None, R_OUT=None, 112 R_FIG=None, L_EXP=None, 113 R_SAVE=None, R_FIGR=None, R_BUF=None, R_BUFR=None, R_BUF_KSH=None, 114 REBUILD_DIR=None, POST_DIR=None, 115 ThreddsPrefix=None, R_GRAF=None, DB=None, 116 IGCM_OUT=None, IGCM_OUT_name=None, rebuild=None, TmpDir=None) : 117 ''' 118 Defines the libIGCM directories - Returns a dictionary or a namespace 119 120 Source : for Spip 121 Possibilities : 122 Local : local (~/Data/IGCMG/...) 123 TGCC_sshfs : TGCC disks mounted via sshfs 124 TGCC_thredds : thredds TGCC via IPSL 125 ('https://thredds-su.ipsl.fr/thredds/dodsC/tgcc_thredds/store/...) 126 127 Exemple of use : 128 libIGCM = libIGCM_sys.config ( .... ) 129 globals().update ( libIGCM ) 130 ''' 131 132 if not Host : Host = Mach (long=False) 133 134 LocalUser = os.environ ['USER'] 135 136 # =========================================================================================== 137 # Reads config.card if available 138 if ConfigCard : 139 # Text existence of ConfigCard 140 if not os.path.exists (ConfigCard ) : 141 raise FileExistsError ( f"File not found : {ConfigCard = }" ) 142 143 ## Creates parser for reading .ini input file 144 MyReader = configparser.ConfigParser (interpolation=configparser.ExtendedInterpolation() ) 145 MyReader.optionxform = str # To keep capitals 146 147 MyReader.read (ConfigCard) 148 149 JobName = MyReader['UserChoices']['JobName'] 150 #----- Short Name of Experiment 151 ExperimentName = MyReader['UserChoices']['ExperimentName'] 152 #----- DEVT TEST PROD 153 SpaceName = MyReader['UserChoices']['SpaceName'] 154 LongName = MyReader['UserChoices']['LongName'] 155 ModelName = MyReader['UserChoices']['ModelName'] 156 TagName = MyReader['UserChoices']['TagName'] 157 158 ### =========================================================================================== 159 ## Part specific to access by OpenDAP/Thredds server 160 161 # =========================================================================================== 162 if Source == 'TGCC_thredds' : 163 if not User and TGCC_User : User = TGCC_User 164 if not Group and TGCC_Group : Group = TGCC_Group 165 166 if not ThreddsPrefix : 167 ThreddsPrefix = 'https://thredds-su.ipsl.fr/thredds/dodsC/tgcc_thredds' 168 169 if not ARCHIVE : ARCHIVE = f'{ThreddsPrefix}/store/{TGCC_User}' 170 if not R_FIG : R_FIG = f'{ThreddsPrefix}/work/{TGCC_User}' 171 if not R_IN : R_IN = f'{ThreddsPrefix}/work/igcmg/IGCM' 172 if not R_GRAF : R_GRAF = f'{ThreddsPrefix}/work/p86mart/GRAF/DATA' 173 174 # =========================================================================================== 175 if Source == 'IDRIS_thredds' : 176 if not User and IDRIS_User : User = IDRIS_User 177 if not Group and IDRIS_Group : Group = IDRIS_Group 178 179 if not ThreddsPrefix : 180 ThreddsPrefix = 'https://thredds-su.ipsl.fr/thredds/catalog/idris_thredds' 181 182 if not ARCHIVE : ARCHIVE = f'{ThreddsPrefix}/store/{IDRIS_User}' 183 if not R_FIG : R_FIG = f'{ThreddsPrefix}/work/{IDRIS_User}' 184 if not R_IN : R_IN = f'{ThreddsPrefix}/work/igcmg/IGCM' 185 if not R_GRAF : R_GRAF = 'https://thredds-su.ipsl.fr/thredds/dodsC/tgcc_thredds/work/p86mart/GRAF/DATA' 186 187 ### =========================================================================================== 188 ## Machine dependant part 189 190 # =========================================================================================== 191 if Host == 'Spip' : 192 if not User : User = LocalUser 193 IGCM_OUT_name = 'IGCM_OUT' 194 if not ARCHIVE : ARCHIVE = os.path.join ( Path.home (), 'Data' ) 195 if not SCRATCHDIR : SCRATCHDIR = os.path.join ( Path.home (), 'Scratch' ) 196 if not R_BUF : R_BUF = os.path.join ( Path.home (), 'Scratch' ) 197 if not R_FIG : R_FIG = os.path.join ( Path.home (), 'Data' ) 198 199 if not STORAGE : STORAGE = ARCHIVE 200 if not R_IN : R_IN = os.path.join ( '/', 'Users', 'marti', 'Data', 'IGCM' ) 201 if not R_GRAF : R_GRAF = os.path.join ( '/', 'Users', 'marti', 'GRAF', 'DATA' ) 202 if not DB : DB = os.path.join ( '/', 'Users', 'marti', 'GRAF', 'DB' ) 203 if not TmpDir : TmpDir = os.path.join ( Path.home (), 'Scratch' ) 204 205 # =========================================================================================== 206 if ( 'Irene' in Host ) or ( 'Rome' in Host ) : 207 208 LocalHome = subprocess.getoutput ( 'ccc_home --ccchome' ) 209 LocalGroup = os.path.basename ( os.path.dirname (LocalHome)) 210 if not User or User == 'marti' : 211 if not TGCC_User : User = LocalUser 212 else : User = TGCC_User 213 214 if not Group : 215 if TGCC_Group : Group = TGCC_Group 216 else : Group = LocalGroup 217 218 IGCM_OUT_name = 'IGCM_OUT' 219 220 if not R_IN : 221 R_IN = os.path.join ( subprocess.getoutput ( 222 'ccc_home --cccwork -d igcmg -u igcmg' ), 'IGCM') 223 if not ARCHIVE : 224 ARCHIVE = subprocess.getoutput ( 225 f'ccc_home --cccstore -u {User} -d {Group}' ) 226 if not STORAGE : 227 STORAGE = subprocess.getoutput ( 228 f'ccc_home --cccwork -u {User} -d {Group}' ) 229 if not SCRATCHDIR : 230 SCRATCHDIR = subprocess.getoutput ( 231 f'ccc_home --cccscratch -u {User} -d {Group}' ) 232 if not R_BUF : 233 R_BUF = subprocess.getoutput ( 234 f'ccc_home --cccscratch -u {User} -d {Group}' ) 235 if not R_FIG : 236 R_FIG = subprocess.getoutput ( 237 f'ccc_home --cccwork -u {User} -d {Group}' ) 238 if not R_GRAF : 239 R_GRAF = os.path.join ( subprocess.getoutput ( 240 'ccc_home --cccwork -d drf -u p86mart'), 'GRAF', 'DATA' ) 241 if not DB : 242 DB = os.path.join ( subprocess.getoutput ( 243 'ccc_home --cccwork -d igcmg -u igcmg'), 'database' ) 244 if not rebuild : 230 245 rebuild = os.path.join ( 231 246 subprocess.getoutput ( 'ccc_home --ccchome -d igcmg -u igcmg' ), 232 247 'Tools', 'irene', 'rebuild_nemo', 'bin', 'rebuild_nemo' ) 233 248 234 if not TmpDir : TmpDir = subprocess.getoutput ( f'ccc_home --cccscratch' ) 235 236 # =========================================================================================== 237 if Host == 'SpiritJ' : 238 if not User : 239 if TGCC_User : User = TGCC_User 240 else : User = LocalUser 241 if not ARCHIVE : 242 ARCHIVE = os.path.join ( '/', 'thredds', 'tgcc', 'store', User ) 243 if not STORAGE : 244 STORAGE = os.path.join ( '/', 'thredds', 'tgcc', 'work' , User ) 245 #if not SCRATCHDIR : SCRATCHDIR = os.path.join ( '/', 'thredds' , 'tgcc', 'store', User ) 246 if not R_IN : 247 R_IN = os.path.join ( '/', 'projsu', 'igcmg', 'IGCM' ) 248 #if not R_GRAF : R_GRAF = os.path.join ('/', 'data', 'omamce', 'GRAF', 'DATA' ) 249 if not R_GRAF : 250 R_GRAF = os.path.join ( '/', 'thredds', 'tgcc', 'work', 'p86mart', 'GRAF', 'DATA' ) 251 if not DB : 252 DB = os.path.join ( '/', 'data', 'igcmg', 'database' ) 253 if not TmpDir : 254 TmpDir = os.path.join ( '/', 'data', LocalUser ) 255 256 # =========================================================================================== 257 if Host == 'SpiritX' : 258 if not User : 259 if TGCC_User : User = TGCC_User 260 else : User = os.environ ['USER'] 261 if not ARCHIVE : 262 ARCHIVE = os.path.join ( '/', 'thredds', 'tgcc', 'store', User ) 263 if not STORAGE : 264 STORAGE = os.path.join ( '/', 'thredds', 'tgcc', 'work' , User ) 265 #if not SCRATCHDIR : 266 # SCRATCHDIR = os.path.join ( '/', 'thredds', 'tgcc', 'store', User ) 267 if not R_IN : 268 R_IN = os.path.join ( '/', 'projsu', 'igcmg', 'IGCM' ) 269 #if not R_GRAF : 270 # R_GRAF = os.path.join ('/', 'data', 'omamce', 'GRAF', 'DATA' ) 271 if not R_GRAF : 272 R_GRAF = os.path.join ( '/', 'thredds', 'tgcc', 'work' , 'p86mart', 'GRAF', 'DATA' ) 273 if not DB : 274 DB = os.path.join ( '/', 'data', 'igcmg', 'database' ) 275 276 # =========================================================================================== 277 if Host == 'Jean-Zay' : 278 if not User : User = os.environ ['USER'] 279 LocalGroup = os.path.basename ( os.path.dirname ( Path.home () )) 280 if not Group : Group = LocalGroup 281 282 if not ARCHIVE : 283 ARCHIVE = os.path.join ( '/', 'gpfsstore' , 'rech', Group, User ) 284 if not STORAGE : 285 STORAGE = os.path.join ( '/', 'gpfswork' , 'rech', Group, User ) 286 if not SCRATCHDIR : 287 SCRATCHDIR = os.path.join ( '/', 'gpfsscratch', 'rech', Group, User ) 288 if not R_FIG : 289 R_FIG = os.path.join ( '/', 'cccwork' , 'rech', Group, User ) 290 if not R_BUF : 291 R_BUF = os.path.join ( '/', 'gpfsscratch', 'rech', Group, User ) 292 if not R_IN : 293 R_IN = os.path.join ( '/', 'gpfswork' , 'rech', 'psl', 'commun', 'IGCM' ) 294 if not R_GRAF : 295 R_GRAF = os.path.join ( '/', 'gpfswork' , 'rech', Group, User, 'GRAF', 'DATA' ) 296 if not DB : 297 DB = os.path.join ( '/', 'gpfswork' , 'rech', 'psl', 'commun', 'database' ) 298 if not rebuild : 299 rebuild = os.path.join ( '/', 'gpfswork', 'rech', 'psl', 'commun', 'Tools', 300 'rebuild', 'modipsl_IOIPSL_PLUS_v2_2_4', 'bin', 'rebuild' ) 301 if not TmpDir : TmpDir = os.path.join ( '/', 'gpfsscratch', 'rech', 302 os.path.basename ( os.path.dirname ( Path.home () )), LocalUser ) 303 304 ### =========================================================================================== 305 ### The construction of the following variables is not machine dependant 306 ### =========================================================================================== 307 if SpaceName == 'TEST' : 308 if SCRATCHDIR and not R_OUT : R_OUT = SCRATCHDIR 309 if SCRATCHDIR and not R_FIG : R_FIG = SCRATCHDIR 310 else : 311 if ARCHIVE and not R_OUT : R_OUT = ARCHIVE 312 if STORAGE and not R_FIG : R_FIG = STORAGE 313 if IGCM_OUT_name : 314 R_OUT = os.path.join ( R_OUT, IGCM_OUT_name ) 315 R_FIG = os.path.join ( R_FIG, IGCM_OUT_name ) 316 317 if SCRATCHDIR and not R_BUF : R_BUF = os.path.join ( SCRATCHDIR, IGCM_OUT_name ) 318 if not IGCM_OUT : IGCM_OUT = R_OUT 319 320 if TagName and SpaceName and ExperimentName and JobName : 321 if not L_EXP : L_EXP = os.path.join ( TagName, SpaceName, ExperimentName, JobName ) 322 323 if R_OUT and not R_SAVE : R_SAVE = os.path.join ( R_OUT , L_EXP ) 249 if not TmpDir : TmpDir = subprocess.getoutput ( f'ccc_home --cccscratch' ) 250 251 # =========================================================================================== 252 if Host == 'SpiritJ' : 253 if not User : 254 if TGCC_User : User = TGCC_User 255 else : User = LocalUser 256 if not ARCHIVE : 257 ARCHIVE = os.path.join ( '/', 'thredds', 'tgcc', 'store', User ) 258 if not STORAGE : 259 STORAGE = os.path.join ( '/', 'thredds', 'tgcc', 'work' , User ) 260 #if not SCRATCHDIR : SCRATCHDIR = os.path.join ( '/', 'thredds' , 'tgcc', 'store', User ) 261 if not R_IN : 262 R_IN = os.path.join ( '/', 'projsu', 'igcmg', 'IGCM' ) 263 #if not R_GRAF : R_GRAF = os.path.join ('/', 'data', 'omamce', 'GRAF', 'DATA' ) 264 if not R_GRAF : 265 R_GRAF = os.path.join ( '/', 'thredds', 'tgcc', 'work', 'p86mart', 'GRAF', 'DATA' ) 266 if not DB : 267 DB = os.path.join ( '/', 'data', 'igcmg', 'database' ) 268 if not TmpDir : 269 TmpDir = os.path.join ( '/', 'data', LocalUser ) 270 271 # =========================================================================================== 272 if Host == 'SpiritX' : 273 if not User : 274 if TGCC_User : User = TGCC_User 275 else : User = os.environ ['USER'] 276 if not ARCHIVE : 277 ARCHIVE = os.path.join ( '/', 'thredds', 'tgcc', 'store', User ) 278 if not STORAGE : 279 STORAGE = os.path.join ( '/', 'thredds', 'tgcc', 'work' , User ) 280 #if not SCRATCHDIR : 281 # SCRATCHDIR = os.path.join ( '/', 'thredds', 'tgcc', 'store', User ) 282 if not R_IN : 283 R_IN = os.path.join ( '/', 'projsu', 'igcmg', 'IGCM' ) 284 if not R_GRAF : 285 R_GRAF = os.path.join ( '/', 'thredds', 'tgcc', 'work' , 'p86mart', 'GRAF', 'DATA' ) 286 if not DB : 287 DB = os.path.join ( '/', 'data', 'igcmg', 'database' ) 288 289 # =========================================================================================== 290 if Host == 'Jean-Zay' : 291 if not User : User = os.environ ['USER'] 292 LocalGroup = os.path.basename ( os.path.dirname ( Path.home () )) 293 if not Group : Group = LocalGroup 294 295 if not ARCHIVE : 296 ARCHIVE = os.path.join ( '/', 'gpfsstore' , 'rech', Group, User ) 297 if not STORAGE : 298 STORAGE = os.path.join ( '/', 'gpfswork' , 'rech', Group, User ) 299 if not SCRATCHDIR : 300 SCRATCHDIR = os.path.join ( '/', 'gpfsscratch', 'rech', Group, User ) 301 if not R_FIG : 302 R_FIG = os.path.join ( '/', 'cccwork' , 'rech', Group, User ) 303 if not R_BUF : 304 R_BUF = os.path.join ( '/', 'gpfsscratch', 'rech', Group, User ) 305 if not R_IN : 306 R_IN = os.path.join ( '/', 'gpfswork' , 'rech', 'psl', 'commun', 'IGCM' ) 307 if not R_GRAF : 308 R_GRAF = os.path.join ( '/', 'gpfswork' , 'rech', Group, User, 'GRAF', 'DATA' ) 309 if not DB : 310 DB = os.path.join ( '/', 'gpfswork' , 'rech', 'psl', 'commun', 'database' ) 311 if not rebuild : 312 rebuild = os.path.join ( '/', 'gpfswork', 'rech', 'psl', 'commun', 'Tools', 313 'rebuild', 'modipsl_IOIPSL_PLUS_v2_2_4', 'bin', 'rebuild' ) 314 if not TmpDir : TmpDir = os.path.join ( '/', 'gpfsscratch', 'rech', 315 os.path.basename ( os.path.dirname ( Path.home () )), LocalUser ) 316 317 ### =========================================================================================== 318 ### The construction of the following variables is not machine dependant 319 ### =========================================================================================== 320 if SpaceName == 'TEST' : 321 if SCRATCHDIR and not R_OUT : R_OUT = SCRATCHDIR 322 if SCRATCHDIR and not R_FIG : R_FIG = SCRATCHDIR 323 else : 324 if ARCHIVE and not R_OUT : R_OUT = ARCHIVE 325 if STORAGE and not R_FIG : R_FIG = STORAGE 324 326 if IGCM_OUT_name : 325 if STORAGE and not R_FIGR : R_FIGR = os.path.join ( STORAGE, IGCM_OUT_name, L_EXP ) 326 else : 327 if STORAGE and not R_FIGR : R_FIGR = os.path.join ( STORAGE, L_EXP ) 328 if R_BUF and not R_BUFR : R_BUFR = os.path.join ( R_BUF , L_EXP ) 329 if R_BUFR and not R_BUF_KSH : R_BUF_KSH = os.path.join ( R_BUFR , 'Out' ) 330 if R_BUF and not REBUILD_DIR : REBUILD_DIR = os.path.join ( R_BUF , L_EXP, 'REBUILD' ) 331 if R_BUF and not POST_DIR : POST_DIR = os.path.join ( R_BUF , L_EXP, 'Out' ) 332 333 ### =========================================================================================== 334 ## Builds the final dictionnary 335 336 libIGCM = { 'TagName' : TagName , 337 'SpaceName' : SpaceName , 338 'ExperimentName': ExperimentName, 339 'JobName' : JobName , 340 'LongName' : LongName , 341 'ModelName' : ModelName , 342 'ShortName' : ShortName, 343 'ConfigCard' : ConfigCard, 344 'ARCHIVE' : ARCHIVE , 345 'STORAGE' : STORAGE , 346 'SCRATCHDIR' : SCRATCHDIR , 347 'R_OUT' : R_OUT , 348 'R_BUF' : R_BUFR , 349 'R_GRAF' : R_GRAF , 350 'DB' : DB , 351 'IGCM_OUT' : IGCM_OUT , 352 'R_SAVE' : R_SAVE , 353 'R_FIGR' : R_FIGR , 354 'R_BUFR' : R_BUFR, 355 'REBUILD_DIR' : REBUILD_DIR, 356 'POST_DIR' : POST_DIR , 357 'R_IN' : R_IN , 358 'Mach' : Host , 359 'Source' : Source , 360 'User' : User , 361 'Group' : Group, 362 'IGCM_OUT_name' : IGCM_OUT_name, 363 'rebuild' : rebuild, 364 'TmpDir' : TmpDir, 365 'TGCC_User' : TGCC_User, 366 'TGCC_Group' : TGCC_Group, 367 'Thredds' : ThreddsPrefix, 368 'IDRIS_User' : IDRIS_User, 369 'IDRIS_Group' : IDRIS_Group } 370 371 if Out in ['namespace', 'space', 'name', 'Space', 'Name', 'Names', 'NameSpace'] : 372 libIGCM = types.SimpleNamespace (**libIGCM) 373 374 return libIGCM 327 R_OUT = os.path.join ( R_OUT, IGCM_OUT_name ) 328 R_FIG = os.path.join ( R_FIG, IGCM_OUT_name ) 329 330 if SCRATCHDIR and not R_BUF : R_BUF = os.path.join ( SCRATCHDIR, IGCM_OUT_name ) 331 if not IGCM_OUT : IGCM_OUT = R_OUT 332 333 if TagName and SpaceName and ExperimentName and JobName : 334 if not L_EXP : L_EXP = os.path.join ( TagName, SpaceName, ExperimentName, JobName ) 335 336 if R_OUT and not R_SAVE : R_SAVE = os.path.join ( R_OUT , L_EXP ) 337 if IGCM_OUT_name : 338 if STORAGE and not R_FIGR : R_FIGR = os.path.join ( STORAGE, IGCM_OUT_name, L_EXP ) 339 else : 340 if STORAGE and not R_FIGR : R_FIGR = os.path.join ( STORAGE, L_EXP ) 341 if R_BUF and not R_BUFR : R_BUFR = os.path.join ( R_BUF , L_EXP ) 342 if R_BUFR and not R_BUF_KSH : R_BUF_KSH = os.path.join ( R_BUFR , 'Out' ) 343 if R_BUF and not REBUILD_DIR : REBUILD_DIR = os.path.join ( R_BUF , L_EXP, 'REBUILD' ) 344 if R_BUF and not POST_DIR : POST_DIR = os.path.join ( R_BUF , L_EXP, 'Out' ) 345 346 ### =========================================================================================== 347 ## Builds the class components 348 349 self.TagName = TagName 350 self.SpaceName = SpaceName 351 self.ExperimentName = ExperimentName 352 self.JobName = JobName 353 self.LongName = LongName 354 self.ModelName = ModelName 355 self.ShortName = ShortName 356 self.ConfigCard = ConfigCard 357 self.ARCHIVE = ARCHIVE 358 self.STORAGE = STORAGE 359 self.SCRATCHDIR = SCRATCHDIR 360 self.R_OUT = R_OUT 361 self.R_BUF = R_BUFR 362 self.R_GRAF = R_GRAF 363 self.DB = DB 364 self.IGCM_OUT = IGCM_OUT 365 self.R_SAVE = R_SAVE 366 self.R_FIGR = R_FIGR 367 self.R_BUFR = R_BUFR 368 self.REBUILD_DIR = REBUILD_DIR 369 self.POST_DIR = POST_DIR 370 self.R_IN = R_IN 371 self.Mach = Host 372 self.Source = Source 373 self.User = User 374 self.Group = Group 375 self.IGCM_OUT_name = IGCM_OUT_name 376 self.rebuild = rebuild 377 self.TmpDir = TmpDir 378 self.TGCC_User = TGCC_User 379 self.TGCC_Group = TGCC_Group 380 self.Thredds = ThreddsPrefix 381 self.IDRIS_User = IDRIS_User 382 self.IDRIS_Group = IDRIS_Group 383 384 -
TOOLS/WATER_BUDGET/lmdz.py
r6665 r6764 178 178 179 179 if verbose : 180 print ( f'{ou_dims =}' 180 print ( f'{ou_dims =}' ) 181 181 print ( f'{new_coords=}' ) 182 182 … … 261 261 return pout 262 262 263 def point2geo (p1d) : 264 '''From 1D (restart type) to 2D 263 def point2geo (p1d, verbose=False, lon=False, lat=False, jpi=0, jpj=0, share_pole=False) : 264 ''' 265 From 1D (..., points_physiques) (restart type) to 2D (..., lat, lon) 265 266 ''' 266 267 math = __mmath__ (p1d) 267 268 268 # Get the dimensions269 # Get the horizontal dimension 269 270 jpn = p1d.shape[-1] 270 271 if len (p1d.shape) > 1 : 272 jpt = p1d.shape[0] 273 else : 274 jpt = 0 275 276 if jpn == 9026 : 277 jpi, jpj = 96, 96 278 if jpn == 17858 : 279 jpi, jpj = 144, 144 280 if jpn == 20306 : 281 jpi, jpj = 144, 143 282 283 if jpt > 0 : 284 p2d = np.zeros ( (jpt, jpj, jpi) ) 285 p2d [:, 1:jpj-1, :] = np.reshape (p1d [:,1:jpn-1], (jpt, jpj-2, jpi) ) 286 p2d [:, 0 , : ] = p1d[:, 0 ] 287 p2d [:, jpj-1, : ] = p1d[:, jpn-1] 288 289 else : 290 p2d = np.zeros ( (jpj, jpi) ) 291 p2d [1:jpj-1, :] = np.reshape (np.float64 (p1d [1:jpn-1]), (jpj-2, jpi) ) 292 p2d [0 , : ] = p1d[ 0 ] 293 p2d [jpj-1, : ] = p1d[jpn-1] 271 # Get other dimension 272 form1 = list(p1d.shape [0:-1]) 273 274 if jpi==0 and jpj>0 : 275 jpi = (jpn-2)//(jpj-2) 276 if jpi>0 and jpj == 0 : 277 jpj = (jpn-2)//jpi +2 278 279 def c_jpn (jpj, jpi) : return jpi*(jpj-2) + 2 280 281 if jpi==0 and jpi==0 : 282 if jpn == c_jpn( 95, 96) : 283 jpj, jpi = 95, 96 284 if jpn == c_jpn (96, 96) : 285 jpj, jpi = 96, 96 286 if jpn == c_jpn (144, 144) : 287 jpj, jpi = 144, 144 288 if jpn == c_jpn (143, 144) : 289 jpj, jpi = 143, 144 290 291 form2 = form1 + [jpj, jpi,] 292 form_shape = form1 + [jpj-2, jpi] 293 p2d = np.empty ( form2 ) 294 295 p2d [..., 1:jpj-1, :] = np.reshape (np.float64 (p1d [..., 1:jpn-1]), form_shape ) 296 if share_pole : 297 for ji in np.arange ( p2d.shape[-1] ) : 298 p2d [..., 0 , ji] = p1d[..., 0 ] / float(jpi) 299 p2d [..., jpj-1 , ji] = p1d[..., jpn-1] / float(jpi) 300 else : 301 for ji in np.arange ( p2d.shape[-1] ) : 302 p2d [..., 0 , ji] = p1d[..., 0 ] 303 p2d [..., jpj-1 , ji] = p1d[..., jpn-1] 294 304 295 305 if math == xr : 296 306 p2d = xr.DataArray (p2d) 297 307 p2d.attrs.update ( p1d.attrs ) 298 p2d = p2d.rename ( {p2d.dims[0]:p1d.dims[0], p2d.dims[-1]:'x', p2d.dims[-2]:'y'} ) 308 for idim in np.arange ( len(p1d.shape [0:-1]) ): 309 dim = p1d.dims[idim] 310 p2d = p2d.rename ( {p2d.dims[idim]:p1d.dims[idim]} ) 311 p2d = p2d.assign_coords ( {p2d.dims[idim]:p1d.coords[dim]} ) 312 p2d = p2d.rename ( {p2d.dims[-1]:'x', p2d.dims[-2]:'y'} ) 313 p2d['x'].attrs.update ( {'axis':'X'} ) 314 p2d['y'].attrs.update ( {'axis':'Y'} ) 315 316 if type(lon) == bool : 317 if lon : 318 lon = np.linspace ( -180, 180, jpi, endpoint=False) 319 #if math == xr : 320 # lon = xr.DataArray ( lon, dims=('lon',), coords=(lon,) ) 321 if type(lat) == bool : 322 if lat : 323 lat = np.linspace ( 90, -90, jpj, endpoint=True) 324 #if math == xr : 325 # lat = xr.DataArray ( lat, dims=('lat',), coords=(lat,) ) 326 if type(lon) != bool and type(lat) != bool : 327 if __mmath__ (lat) == xr : 328 lat_name = lat.name 329 else : 330 lat_name = 'lat' 331 if __mmath__ (lon) == xr : 332 lon_name = lon.name 333 else : 334 lon_name = 'lon' 335 p2d = p2d.rename ( {p2d.dims[-1]:lon_name, p2d.dims[-2]:lat_name} ) 336 p2d = p2d.assign_coords ( {lon_name:lon, lat_name:lat} ) 337 p2d[lon_name].attrs.update ( { 'units':'degrees_east' , 'long_name':'Longitude', 'standard_name':'longitude', 'axis':'X' } ) 338 p2d[lat_name].attrs.update ( { 'units':'degrees_north', 'long_name':'Latitude' , 'standard_name':'latitude' , 'axis':'Y' } ) 299 339 300 340 return p2d 301 341 302 342 def geo2point ( p2d, cumul_poles=False, dim1d='points_physiques' ) : 303 '''From 2D (lat, lon) to 1D (points_phyiques) 343 ''' 344 From 2D (..., lat, lon) to 1D (..., points_phyiques) 304 345 ''' 305 346 math = __mmath__ (p2d) 306 347 # 307 # Get the dimension 308 348 # Get the horizontal dimensions 309 349 (jpj, jpi) = p2d.shape[-2:] 310 311 if len (p2d.shape) > 2 : 312 jpt = p2d.shape[0] 350 jpn = jpi*(jpj-2) + 2 351 352 # Get other dimensions 353 form1 = list(p2d.shape [0:-2]) 354 form2 = form1 + [jpn,] 355 356 p1d = np.empty ( form2 ) 357 form_shape = form1 + [jpn-2,] 358 359 if math == xr : 360 p1d[..., 1:-1] = np.reshape ( np.float64 (p2d[..., 1:-1, :].values).ravel(), form_shape ) 313 361 else : 314 jpt = -1 315 316 jpn = jpi*(jpj-2) + 2 317 318 if jpt == -1 : 319 p1d = np.zeros ( (jpn,) ) 320 p1d[1:-1] = np.float64(p2d[1:-1, :]).ravel() 321 p1d[ 0] = p2d[ 0, 0] 322 p1d[-1] = p2d[-1, 0] 323 324 else : 325 p1d = np.zeros ( (jpt, jpn) ) 326 if math == xr : 327 p1d[:, 1:-1] = np.reshape ( np.float64 (p2d[:, 1:-1, :].values).ravel(), (jpt, jpn-2) ) 328 else : 329 p1d[:, 1:-1] = np.reshape ( np.float64 (p2d[:, 1:-1, :] ).ravel(), (jpt, jpn-2) ) 330 p1d[:, 0 ] = p2d[:, 0, 0] 331 p1d[:, -1 ] = p2d[:, -1, 0] 362 p1d[..., 1:-1] = np.reshape ( np.float64 (p2d[..., 1:-1, :] ).ravel(), form_shape ) 363 p1d[..., 0 ] = p2d[..., 0, 0] 364 p1d[..., -1 ] = p2d[..., -1, 0] 332 365 333 366 if math == xr : 334 367 p1d = xr.DataArray (p1d) 335 368 p1d.attrs.update ( p2d.attrs ) 336 p1d = p1d.rename ( {p1d.dims[0]:p2d.dims[0], p1d.dims[-1]:dim1d} ) 369 if len(p2d.shape [0:-2]) > 0 : 370 for idim in np.arange ( len(p2d.shape [0:-2]) ): 371 dim=p2d.dims[idim] 372 p1d = p1d.rename ( {p1d.dims[idim]:p2d.dims[idim]} ) 373 p1d = p1d.assign_coords ( {p1d.dims[idim] :p2d.coords[dim].values} ) 374 p1d = p1d.rename ( {p1d.dims[-1]:dim1d} ) 337 375 338 376 if cumul_poles : … … 341 379 342 380 return p1d 343 344 def geo3point ( p3d, cumul_poles=False, dim1d='points_physiques' ) :345 '''From 3D (lev, lat, lon) to 2D (lev, points_phyiques)346 '''347 math = __mmath__ (p3d)348 #349 # Get the dimensions350 351 (jpk, jpj, jpi) = p3d.shape[-3:]352 353 if len (p3d.shape) > 3 :354 jpt = p3d.shape[0]355 else :356 jpt = -1357 358 jpn = jpi*(jpj-2) + 2359 360 if jpt == -1 :361 p2d = np.zeros ( (jpk, jpn,) )362 for jk in np.arange (jpk) :363 p2d [jk, :] = geo2point ( p3d [jk,:,:], cumul_poles, dim1d )364 else :365 p2d = np.zeros ( (jpt, jpk, jpn) )366 for jk in np.arange (jpk) :367 p2d [:, jk, :] = geo2point ( p3d [:, jk,:,:], cumul_poles, dim1d )368 369 if math == xr :370 p2d = xr.DataArray (p2d)371 p2d.attrs.update ( p3d.attrs )372 p2d = p2d.rename ( {p2d.dims[-1]:dim1d, p2d.dims[-2]:p3d.dims[-3]} )373 374 return p2d375 381 376 382 def geo2en (pxx, pyy, pzz, glam, gphi) :
Note: See TracChangeset
for help on using the changeset viewer.
