1 | !> \file imposed_regions_rsl.f90 |
---|
2 | !! Module permet d'imposer un scenario des variations locales, par regions - pas par pixel |
---|
3 | !< |
---|
4 | |
---|
5 | !> \namespace imposed_regions_rsl |
---|
6 | !! Module qui permet d'imposer un scenario des variations locales, par regions - pas par pixel |
---|
7 | !! \author A. Quiquet |
---|
8 | !! \date January 2018 |
---|
9 | !! @note used modules |
---|
10 | !! @note - use module3D_phy |
---|
11 | !< |
---|
12 | |
---|
13 | module imposed_regions_rsl |
---|
14 | |
---|
15 | ! 2 fichiers : |
---|
16 | ! i- un fichier texte avec nb_reg colonnes et un certain nombre de lignes, comme fileforc dans forclim |
---|
17 | ! ii- un fichier netcdf qui contient les regions, numerotees en partant de 0. La region 0 n'est pas touchee, niveau eustatique |
---|
18 | |
---|
19 | use module3D_phy,only: nx,ny,time,num_param,num_forc,num_rep_42,dirforcage,dirnameinp,sealevel,sealevel_2d |
---|
20 | use netcdf |
---|
21 | use io_netcdf_grisli |
---|
22 | |
---|
23 | implicit none |
---|
24 | |
---|
25 | |
---|
26 | integer :: nft !< number of lines in the rsl history file |
---|
27 | real,dimension(:),allocatable :: tdate !< time for rsl forcing |
---|
28 | real,dimension(:,:),allocatable :: reg_rsltime !< rsl forcing for a specific region (time,region) |
---|
29 | |
---|
30 | integer nb_reg !< number of regions |
---|
31 | integer,dimension(nx,ny) :: zregions !< mask of the different regions |
---|
32 | |
---|
33 | contains |
---|
34 | |
---|
35 | !> SUBROUTINE: input_rsl |
---|
36 | !! Routine permettant d'introduire des variations locales du niveau marin |
---|
37 | !! en anomalie par rapport au niveau eustatique |
---|
38 | !> |
---|
39 | subroutine input_rsl |
---|
40 | |
---|
41 | implicit none |
---|
42 | |
---|
43 | character(len=8) :: control ! label to check clim. forc. file (filin) is usable |
---|
44 | character(len=150) :: filin ! the file to be read |
---|
45 | character(len=150) :: filinreg ! the file to be read |
---|
46 | character(len=100) :: filforc ! the history of rsl, a row = a region |
---|
47 | character(len=100) :: filreg ! the file of the different regions |
---|
48 | |
---|
49 | |
---|
50 | integer :: err ! recuperation erreur |
---|
51 | |
---|
52 | real*8, dimension(:,:), pointer :: tab !< tableau 2d real pointe |
---|
53 | |
---|
54 | integer i |
---|
55 | |
---|
56 | namelist/imposed_regions_rsl/nb_reg,filreg,filforc |
---|
57 | |
---|
58 | |
---|
59 | rewind(num_param) |
---|
60 | read(num_param,imposed_regions_rsl) |
---|
61 | |
---|
62 | write(num_rep_42,'(A)')'!___________________________________________________________' |
---|
63 | write(num_rep_42,'(A)') '&imposed_regions_rsl ! nom du bloc ' |
---|
64 | write(num_rep_42,*) |
---|
65 | write(num_rep_42,*) 'nb_reg = ', nb_reg |
---|
66 | write(num_rep_42,'(A,A,A)') 'filreg = "',trim(filreg),'"' |
---|
67 | write(num_rep_42,'(A,A,A)') 'filforc = "',trim(filforc),'"' |
---|
68 | write(num_rep_42,*)'/' |
---|
69 | write(num_rep_42,*) |
---|
70 | |
---|
71 | filin=trim(dirforcage)//trim(filforc) |
---|
72 | |
---|
73 | open(num_forc,file=filin,status='old') |
---|
74 | |
---|
75 | read(num_forc,*) control,nft |
---|
76 | |
---|
77 | ! Determination of file size (line nb), allocation of perturbation array |
---|
78 | |
---|
79 | if (control.ne.'nb_lines') then |
---|
80 | write(6,*) filin,'indiquer le nb de ligne en debut de fichier:' |
---|
81 | write(6,*) 'le nb de lignes et le label de control nb_lines' |
---|
82 | stop |
---|
83 | endif |
---|
84 | |
---|
85 | if (.not.allocated(tdate)) then |
---|
86 | allocate(tdate(nft),stat=err) |
---|
87 | if (err/=0) then |
---|
88 | print *,"erreur a l'allocation du tableau Tdate dans imposed_regions_rsl",err |
---|
89 | stop 4 |
---|
90 | end if |
---|
91 | end if |
---|
92 | |
---|
93 | if (.not.allocated(reg_rsltime)) then |
---|
94 | allocate(reg_rsltime(nb_reg,nft),stat=err) |
---|
95 | if (err/=0) then |
---|
96 | print *,"erreur a l'allocation du tableau reg_rsltime dans imposed_regions_rsl",err |
---|
97 | stop 4 |
---|
98 | end if |
---|
99 | end if |
---|
100 | |
---|
101 | do i=1,nft |
---|
102 | read(num_forc,*) tdate(i),reg_rsltime(1:nb_reg,i) |
---|
103 | end do |
---|
104 | close(num_forc) |
---|
105 | |
---|
106 | filinreg=trim(dirnameinp)//trim(filreg) |
---|
107 | call read_ncdf_var('z',filinreg,tab) |
---|
108 | zregions(:,:) = tab(:,:) |
---|
109 | if (minval(zregions(:,:)).ne.0) then |
---|
110 | write(*,*) "Error in imposed_regions_rsl: should be numbered from 0 to nb_reg (min)" |
---|
111 | endif |
---|
112 | if (maxval(zregions(:,:)).ne.nb_reg) then |
---|
113 | write(*,*) "Error in imposed_regions_rsl: should be numbered from 0 to nb_reg (max)" |
---|
114 | endif |
---|
115 | |
---|
116 | end subroutine input_rsl |
---|
117 | |
---|
118 | |
---|
119 | !> SUBROUTINE: rsl |
---|
120 | !! Routine |
---|
121 | !> |
---|
122 | subroutine rsl |
---|
123 | |
---|
124 | integer i,j |
---|
125 | |
---|
126 | real,dimension(0:nb_reg) :: reg_rsl |
---|
127 | |
---|
128 | reg_rsl(0) = sealevel |
---|
129 | |
---|
130 | if(time.lt.tdate(1)) then |
---|
131 | reg_rsl(1:nb_reg)=reg_rsltime(:,1) |
---|
132 | |
---|
133 | else if (time.ge.tdate(nft)) then |
---|
134 | reg_rsl(1:nb_reg)=reg_rsltime(:,nft) |
---|
135 | |
---|
136 | else |
---|
137 | do i=1,nft-1 |
---|
138 | if((time.ge.tdate(i)).and.(time.lt.tdate(i+1))) then ! entre i et i+1 : cas general |
---|
139 | reg_rsl(1:nb_reg)=reg_rsltime(:,i)+(reg_rsltime(:,i+1)-reg_rsltime(:,i))* & |
---|
140 | (time-tdate(i))/(tdate(i+1)-tdate(i)) |
---|
141 | goto 100 |
---|
142 | endif |
---|
143 | end do |
---|
144 | endif |
---|
145 | 100 continue |
---|
146 | |
---|
147 | |
---|
148 | reg_rsl(0) = sealevel |
---|
149 | |
---|
150 | do j=1,ny |
---|
151 | do i=1,nx |
---|
152 | |
---|
153 | sealevel_2d(i,j)=reg_rsl(zregions(i,j)) |
---|
154 | |
---|
155 | enddo |
---|
156 | enddo |
---|
157 | |
---|
158 | |
---|
159 | end subroutine rsl |
---|
160 | |
---|
161 | end module imposed_regions_rsl |
---|