Ignore:
Timestamp:
02/23/16 15:30:13 (8 years ago)
Author:
dumas
Message:

Merge branche iLOVECLIM sur rev 51

Location:
branches/iLoveclim
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • branches/iLoveclim

  • branches/iLoveclim/SOURCES/Hemin40_files/lect-hemin40_mod.f90

    r27 r52  
    22 
    33  use module3D_phy 
     4  use interface_input 
     5  use io_netcdf_grisli 
     6 
     7  implicit none 
    48   
    5     character(len=50) :: FILE1 
    6     character(len=50) :: FILE2 
    7     character(len=80) :: filin 
    8     real, dimension(nx,ny,5) :: bidon          ! pour l'appel a courbure 
     9  character(len=100) :: topo_dep       ! Topo de départ 
     10  character(len=100) :: topo_ref       ! Topo de référence 
     11  character(len=100) :: grid_topo      ! fichier grille 
     12  character(len=100) :: ghf_fich       ! fichier grille 
     13  character(len=80) :: filin 
     14  real, dimension(nx,ny,5) :: bidon          ! pour l'appel a courbure 
     15  character(len=100) :: file_ncdf      !< fichier netcdf issue des fichiers .dat 
    916 
    1017contains 
    … …  
    1219subroutine input_topo 
    1320 
    14 namelist/topo_file/file1,file2 
    15 rewind(num_param)        ! pour revenir au debut du fichier param_list.dat 
    16 read(num_param,topo_file) 
     21  integer :: ios 
     22 
     23  namelist/topo_file/topo_ref,topo_dep,grid_topo,ghf_fich 
     24  rewind(num_param)        ! pour revenir au debut du fichier param_list.dat 
     25  read(num_param,topo_file) 
    1726! formats pour les ecritures dans 42 
    1827428 format(A) 
    19 write(num_rep_42,428)'!___________________________________________________________'  
    20 write(num_rep_42,428) '&topo_file                                  ! nom du bloc ' 
    21 write(num_rep_42,*) 
    22 write(num_rep_42,*) 'file1 = ', file1 
    23 write(num_rep_42,*) 'file2 = ', file2 
    24 write(num_rep_42,*)'/'  
    25 write(num_rep_42,428) '! file1 : topo de depart' 
    26 write(num_rep_42,428) '! file2 : topo de reference'              
    27 write(num_rep_42,*) 
     28  write(num_rep_42,428)'!___________________________________________________________'  
     29  write(num_rep_42,428) '&topo_file                                  ! nom du bloc ' 
     30  write(num_rep_42,*) 
     31  write(num_rep_42,'(A,A)') 'topo_ref = ', topo_ref 
     32  write(num_rep_42,'(A,A)') 'topo_dep = ', topo_dep 
     33  write(num_rep_42,'(A,A)') 'grid_topo =', grid_topo 
     34  write(num_rep_42,'(A,A)') 'ghf_fich = ', ghf_fich 
     35  write(num_rep_42,*)'/'                       
     36  write(num_rep_42,428) '! topo_ref= topo ref isostasie' 
     37  write(num_rep_42,428) '! topo_dep= topo de depart' 
     38  write(num_rep_42,428) '! grid_topo : fichier i,j,x,y,lon,lat' 
     39  write(num_rep_42,428) '! ghf_fich  : fichier flux geothermique'   
     40  write(num_rep_42,*) 
     41 
     42  topo_ref=trim(dirnameinp)//trim(topo_ref)  
     43  topo_dep=trim(dirnameinp)//trim(topo_dep)  
     44  grid_topo=trim(dirnameinp)//trim(grid_topo) 
     45  ghf_fich=trim(dirnameinp)//trim(ghf_fich) 
     46 
     47!write(num_rep_42,*) 'file1 = ', file1 
     48!write(num_rep_42,*) 'file2 = ', file2 
     49!write(num_rep_42,*)'/'  
     50!write(num_rep_42,428) '! file1 : topo de depart' 
     51!write(num_rep_42,428) '! file2 : topo de reference'              
     52!write(num_rep_42,*) 
    2853 
    2954!====================================== La reponse est 42 =========== 
    … …  
    4166 
    4267       
    43 ! lecture adaptee aux fichiers intercomparaison EISMINT 
    44        nxx=nx 
    45        nyy=ny 
     68!!$! lecture adaptee aux fichiers intercomparaison EISMINT 
     69!!$       nxx=nx 
     70!!$       nyy=ny 
    4671 
    47 !     lecture de la topo actuelle 
    48 !     --------------------------- 
    49      open (20,file=TRIM(DIRNAMEINP)//file2,status='old') 
    50         
    51      read(20,'(A80)') TITRE 
    52      read(20,*) NI,NJ,NXX,NYY,STEP 
    53      read(20,*) 
    54          do J=1,ny  
    55           do I=1,nx 
    56              read (20,*)  S0(I,J),H0(I,J),BSOC0(I,J) 
    57              S0(i,j)=max(S0(i,j),0.) 
    58           end do 
    59         end do 
    60      close(20) 
     72!!$!     lecture de la topo actuelle 
     73!!$!     --------------------------- 
     74!!$     open (20,file=TRIM(DIRNAMEINP)//file2,status='old') 
     75!!$        
     76!!$     read(20,'(A80)') TITRE 
     77!!$     read(20,*) NI,NJ,NXX,NYY,STEP 
     78!!$     read(20,*) 
     79!!$         do J=1,ny  
     80!!$          do I=1,nx 
     81!!$             read (20,*)  S0(I,J),H0(I,J),BSOC0(I,J) 
     82!!$             S0(i,j)=max(S0(i,j),0.) 
     83!!$          end do 
     84!!$        end do 
     85!!$     close(20) 
     86!!$ 
     87!!$      
     88!!$!     lecture de la topo de depart 
     89!!$!     --------------------------- 
     90!!$     open (20,file=TRIM(DIRNAMEINP)//file1,status='old') 
     91!!$!         open (20,file='../INPUT-DATA/hemin.g50') 
     92!!$     read(20,'(A80)') TITRE 
     93!!$     read(20,*) NI,NJ,NXX,NYY,STEP 
     94!!$     read(20,*) 
     95!!$         do J=1,ny  
     96!!$          do I=1,nx 
     97!!$             read (20,*) S(I,J),H(I,J),BSOC(I,J) 
     98!!$          end do 
     99!!$        end do 
     100!!$     close(20) 
    61101 
    62       
    63 !     lecture de la topo de depart 
    64 !     --------------------------- 
    65      open (20,file=TRIM(DIRNAMEINP)//file1,status='old') 
    66 !         open (20,file='../INPUT-DATA/hemin.g50') 
    67      read(20,'(A80)') TITRE 
    68      read(20,*) NI,NJ,NXX,NYY,STEP 
    69      read(20,*) 
    70          do J=1,ny  
    71           do I=1,nx 
    72              read (20,*) S(I,J),H(I,J),BSOC(I,J) 
    73           end do 
    74         end do 
    75      close(20) 
     102 
     103!     lecture de la topo de référence 
     104    call lect_input(1,'Bsoc',1,Bsoc0,topo_ref,file_ncdf)    ! socle 
     105    call lect_input(1,'S',1,S0,topo_ref,file_ncdf)          ! surface 
     106    call lect_input(1,'H',1,H0,topo_ref,file_ncdf)          ! epaisseur 
     107 
     108!     lecture de la topo de départ 
     109    call lect_input(1,'Bsoc',1,Bsoc,topo_dep,file_ncdf)    ! socle 
     110    call lect_input(1,'S',1,S,topo_dep,file_ncdf)          ! surface 
     111    call lect_input(1,'H',1,H,topo_dep,file_ncdf)          ! epaisseur 
     112 
    76113 
    77114! calcul des courbures du socle 
    78  
    79115     call courbure(nx,ny,dx,Bsoc,bidon(:,:,1),bidon(:,:,2),bidon(:,:,3), & 
    80116          bidon(:,:,4),socle_cry,bidon(:,:,5)) 
    … …  
    83119! lecture des coordonnées geographiques 
    84120 
    85     filin=TRIM(DIRNAMEINP)//'coord-nord-40km.dat' 
     121!    filin=TRIM(DIRNAMEINP)//grid_topo 
    86122 
    87123! les coordonnees sont calculees en °dec avec GMT, 
    88124! les longitudes sont comprises entre -180 et +180 (negative a l'Ouest de 
    89125! Greenwich et positive a l'Est) 
    90     open(unit=2004,file=filin,iostat=ios) 
     126    open(unit=2004,file=grid_topo,iostat=ios) 
    91127    do k=1,nx*ny 
    92128       read(2004,*) i,j,XCC(i,j),YCC(i,j),XLONG(i,j),YLAT(i,j) 
    … …  
    101137                          
    102138! lecture du flux geothermique de Shapiro 
    103     open(88,file=TRIM(DIRNAMEINP)//'ijphi_hemin40.dat') 
     139!    open(88,file=TRIM(DIRNAMEINP)//'ijphi_hemin40.dat') 
     140!    open(88,file=ghf_fich) 
    104141 
    105     write(42,*) 'flux geothermique Shapiro : ',TRIM(DIRNAMEINP)//'ijphi_hemin40.dat' 
     142!!$    write(42,*) 'flux geothermique Shapiro : ',TRIM(DIRNAMEINP)//'ijphi_hemin40.dat' 
     143!!$    do k=1,nx*ny 
     144!!$       read(88,*) i,j,ghf(i,j) 
     145!!$    end do 
     146!!$    close(88) 
    106147 
    107     do k=1,nx*ny 
    108        read(88,*) i,j,ghf(i,j) 
    109 !        print*, i,j,ghf(i,j) 
    110     end do 
    111     close(88) 
     148    call lect_input(1,'ghf',1,ghf,ghf_fich,file_ncdf) 
     149 
    112150! pour passer les flux des mW/m2 au J/m2/an       
    113151    ghf(:,:)=-SECYEAR/1000.*ghf(:,:) 
Note: See TracChangeset for help on using the changeset viewer.