19 MODULE grid_rect_class
25 DOUBLE PRECISION :: xmin
26 DOUBLE PRECISION :: xmax
27 DOUBLE PRECISION :: ymin
28 DOUBLE PRECISION :: ymax
29 DOUBLE PRECISION :: dx
30 DOUBLE PRECISION :: dy
31 INTEGER :: component_flag
35 MODULE PROCEDURE grid_rect_delete
39 MODULE PROCEDURE grid_rect_get_val
43 MODULE PROCEDURE grid_rect_set_val
47 MODULE PROCEDURE grid_rect_copy
50 INTERFACE OPERATOR(==)
51 MODULE PROCEDURE grid_rect_eq
55 MODULE PROCEDURE grid_rect_write_unit
59 MODULE PROCEDURE grid_rect_read_unit
63 MODULE PROCEDURE grid_rect_display
67 PRIVATE grid_rect_delete, grid_rect_get_val, &
68 grid_rect_set_val, grid_rect_copy, grid_rect_eq, &
69 grid_rect_read_unit, grid_rect_write_unit, grid_rect_display
73 FUNCTION grid_rect_new(xmin, xmax, ymin, ymax, dx, dy, component_flag)
RESULT(this)
74 DOUBLE PRECISION,
INTENT(in),
OPTIONAL :: xmin, xmax, ymin, ymax
75 DOUBLE PRECISION,
INTENT(in),
OPTIONAL :: dx, dy
78 INTEGER,
INTENT(in),
OPTIONAL :: component_flag
80 TYPE(grid_rect) :: this
82 this%xmin = optio_d(xmin)
83 this%ymin = optio_d(ymin)
84 this%xmax = optio_d(xmax)
85 this%ymax = optio_d(ymax)
88 this%component_flag = optio_l(component_flag)
90 END FUNCTION grid_rect_new
93 SUBROUTINE grid_rect_delete(this)
94 TYPE(grid_rect),
INTENT(inout) :: this
102 this%component_flag = imiss
104 END SUBROUTINE grid_rect_delete
107 SUBROUTINE grid_rect_get_val(this, xmin, xmax, ymin, ymax, dx, dy, component_flag)
108 TYPE(grid_rect),
INTENT(in) :: this
109 DOUBLE PRECISION,
INTENT(out),
OPTIONAL :: xmin, xmax, ymin, ymax
110 DOUBLE PRECISION,
INTENT(out),
OPTIONAL :: dx, dy
113 INTEGER,
INTENT(out),
OPTIONAL :: component_flag
115 IF (
PRESENT(xmin))
THEN 118 IF (
PRESENT(ymin))
THEN 121 IF (
PRESENT(xmax))
THEN 124 IF (
PRESENT(ymax))
THEN 127 IF (
PRESENT(dx))
THEN 130 IF (
PRESENT(dy))
THEN 133 IF (
PRESENT(component_flag))
THEN 134 component_flag = this%component_flag
137 END SUBROUTINE grid_rect_get_val
140 SUBROUTINE grid_rect_set_val(this, xmin, xmax, ymin, ymax, &
141 dx, dy, component_flag)
142 TYPE(grid_rect),
INTENT(inout) :: this
143 DOUBLE PRECISION,
INTENT(in),
OPTIONAL :: xmin, xmax, ymin, ymax
144 DOUBLE PRECISION,
INTENT(in),
OPTIONAL :: dx, dy
147 INTEGER,
INTENT(in),
OPTIONAL :: component_flag
150 IF (
PRESENT(xmin))
THEN 153 IF (
PRESENT(ymin))
THEN 156 IF (
PRESENT(xmax))
THEN 159 IF (
PRESENT(ymax))
THEN 162 IF (
PRESENT(dx))
THEN 165 IF (
PRESENT(dy))
THEN 168 IF (
PRESENT(component_flag))
THEN 169 this%component_flag = component_flag
172 END SUBROUTINE grid_rect_set_val
175 SUBROUTINE grid_rect_copy(this, that)
176 TYPE(grid_rect),
INTENT(in) :: this
177 TYPE(grid_rect),
INTENT(out) :: that
181 END SUBROUTINE grid_rect_copy
184 ELEMENTAL FUNCTION grid_rect_eq(this, that)
RESULT(res)
185 TYPE(grid_rect),
INTENT(in) :: this
186 TYPE(grid_rect),
INTENT(in) :: that
191 res = (this%xmin == that%xmin .AND. this%xmax == that%xmax .AND. &
192 this%ymin == that%ymin .AND. this%ymax == that%ymax .AND. &
193 this%dx == that%dx .AND. this%dy == that%dy .AND. &
194 this%component_flag == that%component_flag)
196 END FUNCTION grid_rect_eq
203 SUBROUTINE grid_rect_read_unit(this, unit)
204 TYPE(grid_rect),
INTENT(out) :: this
205 INTEGER,
INTENT(in) :: unit
207 CHARACTER(len=40) :: form
209 INQUIRE(unit, form=form)
210 IF (form ==
'FORMATTED')
THEN 211 READ(unit,*)this%xmin,this%ymin,this%xmax,this%ymax,this%dx,this%dy,this%component_flag
213 READ(unit)this%xmin,this%ymin,this%xmax,this%ymax,this%dx,this%dy,this%component_flag
216 END SUBROUTINE grid_rect_read_unit
223 SUBROUTINE grid_rect_write_unit(this, unit)
224 TYPE(grid_rect),
INTENT(in) :: this
225 INTEGER,
INTENT(in) :: unit
227 CHARACTER(len=40) :: form
229 INQUIRE(unit, form=form)
230 IF (form ==
'FORMATTED')
THEN 231 WRITE(unit,*)this%xmin,this%ymin,this%xmax,this%ymax,this%dx,this%dy,this%component_flag
233 WRITE(unit)this%xmin,this%ymin,this%xmax,this%ymax,this%dx,this%dy,this%component_flag
236 END SUBROUTINE grid_rect_write_unit
240 SUBROUTINE grid_rect_display(this)
241 TYPE(grid_rect),
INTENT(in) :: this
243 print*,
"xFirst",this%xmin
244 print*,
"xLast ",this%xmax
245 print*,
"yFirst",this%ymin
246 print*,
"yLast ",this%ymax
247 print*,
"dx, dy",this%dx,this%dy
248 print*,
"componentFlag",this%component_flag
250 END SUBROUTINE grid_rect_display
256 SUBROUTINE grid_rect_coordinates(this, x, y)
257 TYPE(grid_rect),
INTENT(in) :: this
258 DOUBLE PRECISION,
INTENT(out) :: x(:,:)
259 DOUBLE PRECISION,
INTENT(out) :: y(:,:)
261 DOUBLE PRECISION :: dx, dy
262 INTEGER :: nx, ny, i, j
268 IF (
SIZE(y,1) /= nx .OR.
SIZE(y,2) /= ny)
THEN 275 CALL grid_rect_steps(this, nx, ny, dx, dy)
277 x(:,:) = reshape((/ ((this%xmin+(dx*dble(i)), i=0,nx-1), j=0,ny-1) /),&
279 y(:,:) = reshape((/ ((this%ymin+(dy*dble(j)), i=0,nx-1), j=0,ny-1) /),&
282 END SUBROUTINE grid_rect_coordinates
286 SUBROUTINE grid_rect_steps(this, nx, ny, dx, dy)
287 TYPE(grid_rect),
INTENT(in) :: this
288 INTEGER,
INTENT(in) :: nx
289 INTEGER,
INTENT(in) :: ny
290 DOUBLE PRECISION,
INTENT(out) :: dx
291 DOUBLE PRECISION,
INTENT(out) :: dy
293 IF (
c_e(nx) .AND.
c_e(this%xmax) .AND.
c_e(this%xmin) .AND. &
294 c_e(nx) .AND. nx > 1)
THEN 295 dx = (this%xmax - this%xmin)/dble(nx - 1)
299 IF (
c_e(ny) .AND.
c_e(this%ymax) .AND.
c_e(this%ymin) .AND. &
300 c_e(ny) .AND. ny > 1)
THEN 301 dy = (this%ymax - this%ymin)/dble(ny - 1)
306 END SUBROUTINE grid_rect_steps
310 SUBROUTINE grid_rect_setsteps(this, nx, ny)
311 TYPE(grid_rect),
INTENT(inout) :: this
312 INTEGER,
INTENT(in) :: nx
313 INTEGER,
INTENT(in) :: ny
315 CALL grid_rect_steps(this, nx, ny, this%dx, this%dy)
317 END SUBROUTINE grid_rect_setsteps
319 END MODULE grid_rect_class
Function to check whether a value is missing or not.
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Definitions of constants and functions for working with missing values.