libsim Versione 7.2.1
|
◆ grid_dim_write_unit()
This method writes on a Fortran file unit the contents of the object this. The record can successively be read by the ::read_unit method. The method works both on formatted and unformatted files.
Definizione alla linea 378 del file grid_dim_class.F90. 379! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
380! authors:
381! Davide Cesari <dcesari@arpa.emr.it>
382! Paolo Patruno <ppatruno@arpa.emr.it>
383
384! This program is free software; you can redistribute it and/or
385! modify it under the terms of the GNU General Public License as
386! published by the Free Software Foundation; either version 2 of
387! the License, or (at your option) any later version.
388
389! This program is distributed in the hope that it will be useful,
390! but WITHOUT ANY WARRANTY; without even the implied warranty of
391! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
392! GNU General Public License for more details.
393
394! You should have received a copy of the GNU General Public License
395! along with this program. If not, see <http://www.gnu.org/licenses/>.
396#include "config.h"
397
406IMPLICIT NONE
407
412 INTEGER :: nx
413 INTEGER :: ny
414 DOUBLE PRECISION,POINTER :: lat(:,:)
415 DOUBLE PRECISION,POINTER :: lon(:,:)
417
418INTERFACE delete
419 MODULE PROCEDURE grid_dim_delete
420END INTERFACE
421
422INTERFACE copy
423 MODULE PROCEDURE grid_dim_copy
424END INTERFACE
425
426INTERFACE alloc
427 MODULE PROCEDURE grid_dim_alloc
428END INTERFACE
429
430INTERFACE dealloc
431 MODULE PROCEDURE grid_dim_dealloc
432END INTERFACE
433
434INTERFACE OPERATOR (==)
435 MODULE PROCEDURE grid_dim_eq
436END INTERFACE
437
438INTERFACE write_unit
439 MODULE PROCEDURE grid_dim_write_unit
440END INTERFACE
441
442INTERFACE read_unit
443 MODULE PROCEDURE grid_dim_read_unit
444END INTERFACE
445
446INTERFACE display
447 MODULE PROCEDURE grid_dim_display
448END INTERFACE
449
450PRIVATE grid_dim_delete, grid_dim_copy, grid_dim_alloc, grid_dim_dealloc, &
451 grid_dim_eq, grid_dim_read_unit, grid_dim_write_unit, grid_dim_display
452
453CONTAINS
454
455FUNCTION grid_dim_new(nx, ny) RESULT(this)
456INTEGER, INTENT(in), OPTIONAL :: nx, ny
457
458TYPE(grid_dim) :: this
459
460this%nx = optio_l(nx)
461this%ny = optio_l(ny)
462NULLIFY(this%lat, this%lon)
463
464END FUNCTION grid_dim_new
465
466
467SUBROUTINE grid_dim_delete(this)
468TYPE(grid_dim), INTENT(inout) :: this
469
470CALL dealloc(this)
471this%nx = imiss
472this%ny = imiss
473
474END SUBROUTINE grid_dim_delete
475
476
477SUBROUTINE grid_dim_alloc(this)
478TYPE(grid_dim),INTENT(inout) :: this
479
480IF (ASSOCIATED(this%lon) .AND. ASSOCIATED(this%lat)) THEN
481 IF (SIZE(this%lon, 1) == this%nx .AND. SIZE(this%lon, 2) == this%ny .AND. &
482 SIZE(this%lat, 1) == this%nx .AND. SIZE(this%lat, 2) == this%ny) RETURN
483ENDIF
484CALL dealloc(this)
486 ALLOCATE(this%lon(this%nx, this%ny), this%lat(this%nx, this%ny))
487ENDIF
488
489END SUBROUTINE grid_dim_alloc
490
491
492SUBROUTINE grid_dim_dealloc(this)
493TYPE(grid_dim),INTENT(inout) :: this
494
495IF (ASSOCIATED(this%lon)) DEALLOCATE(this%lon)
496IF (ASSOCIATED(this%lat)) DEALLOCATE(this%lat)
497
498END SUBROUTINE grid_dim_dealloc
499
500
501SUBROUTINE grid_dim_copy(this, that)
502TYPE(grid_dim),INTENT(in) :: this
503TYPE(grid_dim),INTENT(out) :: that
504
505that = grid_dim_new(this%nx, this%ny)
506
507IF (ASSOCIATED(this%lon) .AND. ASSOCIATED(this%lat))THEN
508 CALL alloc(that)
509
510#ifdef DEBUG
511 IF (SIZE(this%lon,1) /= this%nx .OR. SIZE(this%lon,2) /= this%ny) THEN
512 CALL raise_error('grid_dim_copy, dimensioni non valide: '// &
515 ENDIF
516 IF (SIZE(this%lat,1) /= this%nx .OR. SIZE(this%lat,2) /= this%ny) THEN
517 CALL raise_error('grid_dim_copy, dimensioni non valide: '// &
520 ENDIF
521#endif
522
523 that%lon(:,:) = this%lon(:,:)
524 that%lat(:,:) = this%lat(:,:)
525ENDIF
526
527END SUBROUTINE grid_dim_copy
528
529
530ELEMENTAL FUNCTION grid_dim_eq(this, that) RESULT(res)
531TYPE(grid_dim),INTENT(IN) :: this, that
532LOGICAL :: res
533
534res = this%nx == that%nx .and. &
535 this%ny == that%ny
536
537END FUNCTION grid_dim_eq
538
539
544SUBROUTINE grid_dim_read_unit(this, unit)
545TYPE(grid_dim),INTENT(out) :: this
546INTEGER, INTENT(in) :: unit
547
548CHARACTER(len=40) :: form
549LOGICAL :: is_all
550
551INQUIRE(unit, form=form)
552IF (form == 'FORMATTED') THEN
553 READ(unit,*)this%nx,this%ny
554 READ(unit,*)is_all
555 IF (is_all) THEN
556 CALL alloc(this)
557 READ(unit,*)this%lon,this%lat
558 ELSE
559 READ(unit,*)
560 ENDIF
561ELSE
562 READ(unit)this%nx,this%ny
563 READ(unit)is_all
564 IF (is_all) THEN
565 CALL alloc(this)
566 READ(unit)this%lon,this%lat
567 ELSE
568 READ(unit)
569 ENDIF
570ENDIF
571
572END SUBROUTINE grid_dim_read_unit
573
574
579SUBROUTINE grid_dim_write_unit(this, unit)
580TYPE(grid_dim),INTENT(in) :: this
581INTEGER, INTENT(in) :: unit
582
583CHARACTER(len=40) :: form
584LOGICAL :: is_all
585
586INQUIRE(unit, form=form)
587IF (form == 'FORMATTED') THEN
588 WRITE(unit,*)this%nx,this%ny
589 is_all = (ASSOCIATED(this%lon) .AND. ASSOCIATED(this%lat))
590 WRITE(unit,*)is_all
591 IF (is_all) THEN
592 WRITE(unit,*)this%lon,this%lat
593 ELSE
594 WRITE(unit,*)
595 ENDIF
596ELSE
597 WRITE(unit)this%nx,this%ny
598 is_all = (ASSOCIATED(this%lon) .AND. ASSOCIATED(this%lat))
599 WRITE(unit)is_all
600 IF (is_all) THEN
601 WRITE(unit)this%lon,this%lat
602 ELSE
603 WRITE(unit)
604 ENDIF
605ENDIF
606
607END SUBROUTINE grid_dim_write_unit
608
609
611SUBROUTINE grid_dim_display(this)
612TYPE(grid_dim),INTENT(in) :: this
613
614print*,'Number of points along x direction',this%nx
615print*,'Number of points along y direction',this%ny
616
617END SUBROUTINE grid_dim_display
618
Set of functions that return a CHARACTER representation of the input variable. Definition char_utilities.F90:253 Function to check whether a value is missing or not. Definition missing_values.f90:72 Module for defining the extension and coordinates of a rectangular georeferenced grid. Definition grid_dim_class.F90:211 Definitions of constants and functions for working with missing values. Definition missing_values.f90:50 Module for quickly interpreting the OPTIONAL parameters passed to a subprogram. Definition optional_values.f90:28 Derived type describing the extension of a grid and the geographical coordinates of each point. Definition grid_dim_class.F90:221 |