Line data Source code
1 : !-----------------------------------------------------------------------
2 : !$Id$
3 : !===============================================================================
4 : module input_reader
5 :
6 : ! Description:
7 : ! This module is respondsible for the procedures and structures necessary to
8 : ! read in "SAM-Like" case specific files. Currently only the
9 : ! <casename>_sounding.in file is formatted to be used by this module.
10 : !
11 : ! References:
12 : ! None
13 : !---------------------------------------------------------------------------------------------------
14 :
15 : use clubb_precision, only: &
16 : core_rknd ! Variable(s)
17 :
18 : implicit none
19 :
20 : private
21 :
22 : public :: one_dim_read_var, &
23 : read_one_dim_file, &
24 : two_dim_read_var, &
25 : read_two_dim_file, &
26 : fill_blanks_one_dim_vars, &
27 : fill_blanks_two_dim_vars, &
28 : deallocate_one_dim_vars, &
29 : deallocate_two_dim_vars, &
30 : read_x_table, &
31 : read_x_profile, &
32 : get_target_index, &
33 : count_columns
34 :
35 : ! Derived type for representing a rank 1 variable that has been read in by one
36 : ! of the procedures.
37 : type one_dim_read_var
38 :
39 : character(len=30) :: name ! Name of the variable
40 :
41 : character(len=30) :: dim_name ! Name of the dimension that the
42 : ! variable varies along
43 :
44 : real( kind = core_rknd ), dimension(:), allocatable :: values ! Values of that variable
45 :
46 : end type one_dim_read_var
47 :
48 : ! Derived type for representing a rank 2 variable that has been read in by one
49 : ! of the procedures.
50 : type two_dim_read_var
51 :
52 : character(len=30) :: name ! Name of the variable
53 :
54 : character(len=30) :: dim1_name ! Name of one of the dimensions
55 : ! that the variable varies along
56 :
57 : character(len=30) :: dim2_name ! Name of the other variable that
58 : ! the variable varies along
59 :
60 : real( kind = core_rknd ), dimension(:,:), allocatable :: values ! Values of that variable
61 :
62 : end type two_dim_read_var
63 :
64 :
65 : ! Constant Parameter(s)
66 : real( kind = core_rknd ), parameter, private :: &
67 : blank_value = -999.9_core_rknd ! Used to denote if a value is missing from the file
68 :
69 : contains
70 :
71 : !-------------------------------------------------------------------------------------------------
72 0 : subroutine read_two_dim_file( iunit, nCol, filename, read_vars, other_dim )
73 : !
74 : ! Description: This subroutine reads from a file containing data that varies
75 : ! in two dimensions. These are dimensions are typically height
76 : ! and time.
77 : !
78 : !-----------------------------------------------------------------------------------------------
79 : use constants_clubb, only: &
80 : fstderr ! Constant(s)
81 :
82 : use input_names, only: &
83 : time_name ! Constant(s)
84 :
85 : use clubb_precision, only: &
86 : core_rknd ! Variable(s)
87 :
88 : implicit none
89 :
90 : ! External
91 : intrinsic :: trim, index
92 :
93 : ! Input Variable(s)
94 :
95 : integer, intent(in) :: iunit ! File I/O unit
96 :
97 : integer, intent(in) :: nCol ! Number of columns expected in the data file
98 :
99 :
100 : character(len=*), intent(in) :: filename ! Name of the file being read from
101 :
102 : ! Output Variable(s)
103 : type (two_dim_read_var), dimension(nCol),intent(out) :: read_vars ! Structured information
104 : ! from the file
105 :
106 : type (one_dim_read_var), intent(out) :: other_dim ! Structured information
107 : ! on the dimesion not stored in read_vars
108 :
109 : ! Local Variables
110 0 : character(len=30),dimension(nCol) :: names ! Names of variables
111 :
112 : integer nRowI ! Inner row
113 :
114 : integer nRowO ! Outer row
115 :
116 : integer :: k, j, i
117 :
118 : logical :: isComment
119 :
120 : character(len=200) :: tmpline
121 :
122 0 : real( kind = core_rknd ), dimension(nCol) :: tmp
123 :
124 : integer :: input_status ! The status of a read statement
125 :
126 : ! ---- Begin Code ----
127 :
128 : ! First run through, take names and determine how large the data file is.
129 0 : open(unit=iunit, file=trim( filename ), status = 'old', action='read' )
130 :
131 0 : isComment = .true.
132 :
133 : ! Skip all the comments at the top of the file
134 : do while ( isComment )
135 0 : read(iunit,fmt='(A)') tmpline
136 0 : k = index( tmpline, "!" )
137 0 : isComment = .false.
138 0 : if ( k > 0 ) then
139 : isComment = .true.
140 : end if
141 : end do
142 :
143 : ! Go back to the line that wasn't a comment.
144 0 : backspace(iunit)
145 :
146 0 : read(iunit, fmt=*) names
147 :
148 0 : nRowO = 0
149 0 : do while(.true.)
150 0 : read(iunit, *, iostat=input_status) tmp(1), nRowI
151 :
152 : ! If input_status shows an end of data, then exit the loop
153 0 : if( input_status < 0 ) then
154 : exit
155 0 : else if ( input_status > 0 ) then
156 0 : write(fstderr,*) "Error reading data from file: " //trim( filename )
157 0 : error stop "Fatal error input_reader"
158 : end if
159 :
160 0 : if( nRowI < 1 ) then
161 : error stop &
162 0 : "Number of elements must be an integer and greater than zero in two-dim input file."
163 : end if
164 :
165 0 : do k =1, nRowI
166 0 : read(iunit, *) tmp
167 : end do
168 0 : nRowO = nRowO + 1
169 : end do
170 :
171 0 : do i=1, nRowO
172 :
173 0 : backspace(iunit)
174 :
175 0 : do j=1, nRowI
176 :
177 0 : backspace(iunit)
178 :
179 : end do
180 :
181 : end do
182 :
183 0 : backspace(iunit)
184 :
185 : ! Store the names into the structure and allocate accordingly
186 0 : do k =1, nCol
187 0 : read_vars(k)%name = names(k)
188 0 : read_vars(k)%dim1_name = time_name
189 0 : read_vars(k)%dim2_name = names(1)
190 :
191 0 : allocate( read_vars(k)%values(nRowI, nRowO) )
192 : end do
193 :
194 0 : other_dim%name = time_name
195 0 : other_dim%dim_name = time_name
196 :
197 0 : allocate( other_dim%values(nRowO) )
198 :
199 : ! Read in the data again to the newly allocated arrays
200 0 : do k=1, nRowO
201 0 : read(iunit,*) other_dim%values(k)
202 0 : do j=1, nRowI
203 0 : read(iunit,*) ( read_vars(i)%values(j,k), i=1, nCol)
204 : end do
205 : end do
206 :
207 0 : close(iunit)
208 :
209 : ! Eliminate a compiler warning
210 : if ( .false. ) print *, tmp
211 :
212 0 : return
213 : end subroutine read_two_dim_file
214 :
215 : !------------------------------------------------------------------------------------------------
216 0 : subroutine read_one_dim_file( iunit, nCol, filename, &
217 0 : read_vars )
218 : !
219 : ! Description:
220 : ! This subroutine reads from a file containing data that varies
221 : ! in one dimension. The dimension is typically time.
222 : !
223 : ! References:
224 : ! None
225 : !----------------------------------------------------------------------------------------------
226 :
227 : use clubb_precision, only: &
228 : core_rknd ! Variable(s)
229 :
230 : implicit none
231 :
232 : ! External
233 :
234 : intrinsic :: trim, index
235 :
236 : ! Input Variable(s)
237 :
238 : integer, intent(in) :: iunit ! I/O unit
239 :
240 : integer, intent(in) :: nCol ! Number of columns expected in the data file
241 :
242 : character(len=*), intent(in) :: filename ! Name of the file being read from
243 :
244 : ! Output Variable(s)
245 :
246 : type (one_dim_read_var), dimension(nCol),intent(out) :: &
247 : read_vars ! Structured information from the file
248 :
249 : ! Local Variable(s)
250 0 : character(len=30),dimension(nCol) :: names
251 :
252 : character(len=200) :: tmpline
253 :
254 : integer nRow
255 :
256 : integer :: k, j
257 :
258 0 : real( kind = core_rknd ), dimension(nCol) :: tmp
259 :
260 : logical :: isComment
261 :
262 : integer :: input_status ! The status of a read statement
263 :
264 : ! ---- Begin Code ----
265 :
266 0 : isComment = .true.
267 :
268 : ! First run through, take names and determine how large the data file is.
269 0 : open(unit=iunit, file=trim( filename ), status = 'old' )
270 :
271 : ! Skip all the comments at the top of the file
272 : do while(isComment)
273 0 : read(iunit,fmt='(A)') tmpline
274 0 : k = index( tmpline, "!" )
275 0 : isComment = .false.
276 0 : if(k > 0) then
277 : isComment = .true.
278 : end if
279 : end do
280 :
281 : ! Go back to the line that wasn't a comment.
282 0 : backspace(iunit)
283 :
284 0 : read(iunit, fmt=*) names
285 :
286 : ! Count up the number of rows
287 0 : nRow = 0
288 0 : do while(.true.)
289 0 : read(iunit, *, iostat=input_status) tmp
290 :
291 : ! If input_status shows an end of file, exit the loop
292 0 : if( input_status < 0 ) then
293 : exit
294 : end if
295 :
296 0 : nRow = nRow+1
297 : end do
298 :
299 : ! Rewind that many rows
300 0 : do k = 0, nRow
301 0 : backspace(iunit)
302 : end do
303 :
304 : ! Store the names into the structure and allocate accordingly
305 0 : do k = 1, nCol
306 0 : read_vars(k)%name = names(k)
307 0 : read_vars(k)%dim_name = names(1)
308 0 : allocate( read_vars(k)%values(nRow) )
309 : end do
310 :
311 : ! Read in the data again to the newly allocated arrays
312 0 : do k=1, nRow
313 0 : read(iunit,*) ( read_vars(j)%values(k), j=1, nCol)
314 : end do
315 :
316 0 : close(iunit)
317 :
318 : ! Avoiding compiler warning
319 : if ( .false. ) print *, tmp
320 :
321 0 : return
322 :
323 : end subroutine read_one_dim_file
324 :
325 : !------------------------------------------------------------------------------------------------
326 0 : subroutine fill_blanks_one_dim_vars( num_vars, one_dim_vars )
327 : !
328 : ! Description:
329 : ! This subroutine fills in the blank spots (signified by constant blank_value)
330 : ! with values linearly interpolated using the first element of the array as a
331 : ! guide.
332 : !
333 : ! References:
334 : ! None
335 : !----------------------------------------------------------------------------------------------
336 :
337 : use clubb_precision, only: &
338 : core_rknd ! Variable(s)
339 :
340 : implicit none
341 :
342 : ! External
343 : intrinsic :: size
344 :
345 : ! Input Variable(s)
346 : integer, intent(in) :: num_vars ! Number of elements in one_dim_vars
347 :
348 : ! Input/Output Variable(s)
349 : type(one_dim_read_var), dimension(num_vars), intent(inout) :: &
350 : one_dim_vars ! Read data that may have gaps.
351 :
352 : ! Local variable(s)
353 : integer :: i
354 :
355 : ! ---- Begin Code ----
356 :
357 0 : do i=1, num_vars
358 0 : one_dim_vars(i)%values = linear_fill_blanks( size( one_dim_vars(i)%values ), &
359 0 : one_dim_vars(1)%values, one_dim_vars(i)%values, &
360 0 : 0.0_core_rknd )
361 : end do
362 :
363 0 : return
364 :
365 : end subroutine fill_blanks_one_dim_vars
366 :
367 : !------------------------------------------------------------------------------------------------
368 0 : subroutine fill_blanks_two_dim_vars( num_vars, other_dim, two_dim_vars )
369 : !
370 : ! Description:
371 : ! This subroutine fills in the blank spots (signified by the
372 : ! constant blank_value with values linearly interpolated using the first
373 : ! element of the array and the values in the other_dim argument as a guide.
374 : !
375 : ! This is a two step process. First we assume that the other_dim values
376 : ! have no holes, but there are blanks for that variable across that
377 : ! dimension. Then we fill holes across the dimension whose values are first
378 : ! in the array of two_dim_vars.
379 : !
380 : ! Ex. Time is the 'other_dim' and Height in meters is the first element in
381 : ! two_dim_vars.
382 : !
383 : ! References:
384 : ! None
385 : !----------------------------------------------------------------------------------------------
386 :
387 : implicit none
388 :
389 : ! External
390 : intrinsic :: size
391 :
392 : ! Input Variable(s)
393 : integer, intent(in) :: num_vars ! Number of elements in one_dim_vars
394 :
395 : ! Input/Output Variable(s)
396 : type(one_dim_read_var), intent(in) :: other_dim ! Read data
397 :
398 : type(two_dim_read_var), dimension(num_vars), intent(inout) :: &
399 : two_dim_vars ! Read data that may have gaps.
400 :
401 : ! Local variables
402 : integer :: i,j ! Loop iterators
403 :
404 : integer :: &
405 : dim_size, & ! 1st dimension size
406 : other_dim_size ! 2nd dimension size
407 :
408 : ! ---- Begin Code ----
409 :
410 0 : dim_size = size( two_dim_vars(1)%values, 1 )
411 :
412 0 : other_dim_size = size( other_dim%values )
413 :
414 0 : do i=2, num_vars
415 : ! Interpolate along main dim
416 0 : do j=1, other_dim_size
417 0 : two_dim_vars(i)%values(:,j) = linear_fill_blanks( dim_size, &
418 : two_dim_vars(1)%values(:,j), &
419 0 : two_dim_vars(i)%values(:,j), blank_value )
420 : end do ! j = 1 .. other_dim_size
421 :
422 : ! Interpolate along other dim
423 0 : do j=1, dim_size
424 0 : two_dim_vars(i)%values(j,:) = linear_fill_blanks( other_dim_size, &
425 : other_dim%values, &
426 0 : two_dim_vars(i)%values(j,:), blank_value )
427 : end do ! j = 1 .. dim_size
428 :
429 : end do ! i = 2 .. num_vars
430 :
431 0 : return
432 :
433 : end subroutine fill_blanks_two_dim_vars
434 :
435 :
436 : !------------------------------------------------------------------------------------------------
437 0 : function linear_fill_blanks( dim_grid, grid, var, default_value ) &
438 0 : result( var_out )
439 : !
440 : ! Description:
441 : ! This function fills blanks in array var using the grid
442 : ! as a guide. Blank values in var are signified by being
443 : ! less than or equal to the constant blank_value.
444 : !
445 : ! References:
446 : ! None
447 : !-----------------------------------------------------------------------------------------------
448 :
449 : use interpolation, only: zlinterp_fnc
450 :
451 : use clubb_precision, only: &
452 : core_rknd ! Variable(s)
453 :
454 : implicit none
455 :
456 : ! Input Variable(s)
457 : integer, intent(in) :: dim_grid ! Size of grid
458 :
459 : real( kind = core_rknd ), dimension(dim_grid), intent(in) :: &
460 : grid ! Array that var is being interpolated to.
461 :
462 : real( kind = core_rknd ), dimension(dim_grid), intent(in) :: &
463 : var ! Array that may contain gaps.
464 :
465 : real( kind = core_rknd ), intent(in) :: &
466 : default_value ! Default value if entire profile == blank_value
467 :
468 : ! Output Variable(s)
469 : real( kind = core_rknd ), dimension(dim_grid) :: &
470 : var_out ! Return variable
471 :
472 : ! Local Variables
473 0 : real( kind = core_rknd ), dimension(dim_grid) :: temp_grid
474 0 : real( kind = core_rknd ), dimension(dim_grid) :: temp_var
475 :
476 : integer :: i
477 : integer :: amt
478 :
479 : logical :: reversed
480 :
481 : ! ---- Begin Code ----
482 :
483 0 : reversed = .false.
484 :
485 : ! Essentially this code leverages the previously written zlinterp function.
486 : ! A smaller temporary grid and var variable are being created to pass to
487 : ! zlinterp. zlinterp then performs the work of taking the temporary var
488 : ! array and interpolating it to the actual grid array.
489 :
490 0 : amt = 0
491 0 : do i=1, dim_grid
492 0 : if ( var(i) > blank_value ) then
493 0 : amt = amt + 1
494 0 : temp_var(amt) = var(i)
495 0 : temp_grid(amt) = grid(i)
496 : end if
497 0 : if ( i > 1 ) then
498 0 : if ( grid(i) < grid(i-1) ) then
499 0 : reversed = .true.
500 : end if
501 : end if
502 : end do
503 :
504 :
505 0 : if ( amt == 0 ) then
506 0 : var_out = default_value
507 0 : else if (amt < dim_grid) then
508 0 : if ( reversed ) then
509 0 : var_out = zlinterp_fnc( dim_grid, amt, -grid, -temp_grid(1:amt), temp_var(1:amt) )
510 : else
511 0 : var_out = zlinterp_fnc( dim_grid, amt, grid, temp_grid(1:amt), temp_var(1:amt) )
512 : end if
513 : else
514 0 : var_out = var
515 : end if
516 :
517 0 : return
518 0 : end function linear_fill_blanks
519 : !----------------------------------------------------------------------------
520 0 : subroutine deallocate_one_dim_vars( num_vars, &
521 0 : one_dim_vars )
522 : !
523 : ! Description:
524 : ! This subroutine deallocates the pointer stored in
525 : ! one_dim_vars%value for the whole array.
526 : !
527 : !------------------------------------------------------------------------------
528 : implicit none
529 :
530 : ! External functions
531 : intrinsic :: allocated
532 :
533 : ! Input Variable(s)
534 : integer, intent(in) :: num_vars ! Number of elements in one_dim_vars
535 :
536 : type(one_dim_read_var), dimension(num_vars), intent(inout) :: &
537 : one_dim_vars ! Read data that may have gaps.
538 :
539 : ! Local Variable(s)
540 : integer :: i
541 :
542 : ! Begin Code
543 :
544 0 : do i=1, num_vars
545 :
546 0 : if ( allocated( one_dim_vars(i)%values ) ) then
547 :
548 0 : deallocate( one_dim_vars(i)%values )
549 :
550 : end if
551 :
552 : end do ! 1 .. num_vars
553 :
554 0 : return
555 : end subroutine deallocate_one_dim_vars
556 :
557 : !------------------------------------------------------------------------------------------------
558 0 : subroutine deallocate_two_dim_vars( num_vars, two_dim_vars, other_dim )
559 : !
560 : ! Description:
561 : ! This subroutine deallocates the pointer stored in
562 : ! two_dim_vars%value for the whole array
563 : !
564 : ! References:
565 : ! None
566 : !----------------------------------------------------------------------------------------------
567 : implicit none
568 :
569 : ! External Functions
570 : intrinsic :: allocated
571 :
572 : ! Input Variable(s)
573 : integer, intent(in) :: num_vars ! Number of elements in one_dim_vars
574 :
575 : ! Input/Output Variables
576 : type(one_dim_read_var), intent(inout) :: other_dim
577 :
578 : type(two_dim_read_var), dimension(num_vars), intent(inout) :: &
579 : two_dim_vars ! Read data that may have gaps.
580 :
581 : ! Local Variable(s)
582 : integer :: i
583 :
584 : ! ---- Begin Code ----
585 :
586 0 : do i=1, num_vars
587 :
588 0 : if ( allocated( two_dim_vars(i)%values ) ) then
589 :
590 0 : deallocate(two_dim_vars(i)%values)
591 :
592 : end if
593 :
594 : end do
595 :
596 0 : if ( allocated( other_dim%values ) ) then
597 :
598 0 : deallocate(other_dim%values)
599 :
600 : end if
601 :
602 0 : return
603 : end subroutine deallocate_two_dim_vars
604 : !------------------------------------------------------------------------------------------------
605 0 : function read_x_table( nvar, xdim, ydim, target_name, retVars ) result( x )
606 : !
607 : ! Description:
608 : ! Searches for the variable specified by target_name in the
609 : ! collection of retVars. If the function finds the variable then it returns
610 : ! it. If it does not the program using this function will exit gracefully
611 : ! with a warning message.
612 : !
613 : ! References:
614 : ! None
615 : !-----------------------------------------------------------------------------------------------
616 :
617 : use clubb_precision, only: &
618 : core_rknd ! Variable(s)
619 :
620 : use constants_clubb, only: &
621 : fstderr ! Constant(s)
622 :
623 : implicit none
624 :
625 : ! External Functions
626 : intrinsic :: trim
627 :
628 : ! Input Variable(s)
629 : integer, intent(in) :: nvar ! Number of variables in retVars
630 :
631 : integer, intent(in) :: xdim, ydim
632 :
633 : character(len=*), intent(in) :: &
634 : target_name ! Name of the variable that is being searched for
635 :
636 : type(two_dim_read_var), dimension(nvar), intent(in) :: &
637 : retVars ! Collection of data being searched through
638 :
639 : ! Output Variable(s)
640 : real( kind = core_rknd ), dimension( xdim, ydim ) :: x
641 :
642 : ! Local Variables
643 : integer :: i ! Loop iterator
644 :
645 : logical :: l_found
646 :
647 : ! ---- Begin Code ----
648 :
649 0 : l_found = .false.
650 :
651 0 : i = 1
652 :
653 0 : do while( i <= nvar .and. .not. l_found)
654 :
655 0 : if( retVars(i)%name == target_name ) then
656 :
657 0 : l_found = .true.
658 :
659 0 : x = retVars(i)%values
660 :
661 : end if
662 :
663 0 : i=i+1
664 :
665 : end do ! i <= nvar .and. not l_found
666 :
667 0 : if ( .not. l_found ) then
668 :
669 0 : write(fstderr,*) trim( target_name )//" could not be found."
670 :
671 0 : error stop "Fatal error in function read_x_table"
672 :
673 : end if
674 :
675 0 : return
676 :
677 0 : end function read_x_table
678 :
679 :
680 : !------------------------------------------------------------------------------------------------
681 0 : function read_x_profile( nvar, dim_size, target_name, retVars, &
682 0 : input_file ) result( x )
683 : !
684 : ! Description:
685 : ! Searches for the variable specified by target_name in the
686 : ! collection of retVars. If the function finds the variable then it returns
687 : ! it. If it does not the program using this function will exit gracefully
688 : ! with a warning message.
689 : !
690 : ! Modified by Cavyn, June 2010
691 : !----------------------------------------------------------------------------------------------
692 :
693 : use constants_clubb, only: &
694 : fstderr ! Variable for writing to error stream
695 :
696 : use clubb_precision, only: &
697 : core_rknd ! Variable(s)
698 :
699 : implicit none
700 :
701 : ! External Functions
702 : intrinsic :: present, size, trim
703 :
704 : ! Input Variable(s)
705 : integer, intent(in) :: &
706 : nvar, & ! Number of variables in retVars
707 : dim_size ! Size of the array returned
708 :
709 : character(len=*), intent(in) :: &
710 : target_name ! Name of the variable that is being searched for
711 :
712 : type(one_dim_read_var), dimension(nvar), intent(in) :: &
713 : retVars ! Collection being searched
714 :
715 : character(len=*), optional, intent(in) :: &
716 : input_file ! Name of the input file containing the variables
717 :
718 : ! Output Variable(s)
719 : real( kind = core_rknd ), dimension(dim_size) :: x
720 :
721 : ! Local Variables
722 : integer :: i
723 :
724 : ! ---- Begin Code ----
725 :
726 0 : i = get_target_index( nvar, target_name, retVars )
727 :
728 0 : if ( i > 0 ) then
729 0 : x(1:size(retVars(i)%values)) = retVars(i)%values
730 :
731 : else
732 0 : if( present( input_file ) ) then
733 0 : write(fstderr,*) trim( target_name ), ' could not be found. Check the file ', input_file
734 : else
735 0 : write(fstderr,*) trim( target_name ), ' could not be found. Check your sounding.in file.'
736 : end if ! present( input_file )
737 0 : error stop "Fatal error in read_x_profile"
738 :
739 : end if ! target_exists_in_array
740 :
741 0 : return
742 :
743 0 : end function read_x_profile
744 :
745 : !------------------------------------------------------------------------------
746 0 : function get_target_index( nvar, target_name, retVars) result( i )
747 : !
748 : ! Description:
749 : ! Returns the index of the variable specified by target_name in the
750 : ! collection of retVars. Returns -1 if variable does not exist in retVars
751 : !
752 : ! References:
753 : ! None
754 : !
755 : ! Created by Cavyn, July 2010
756 : !----------------------------------------------------------------------------------------------
757 :
758 : implicit none
759 :
760 : ! Input Variable(s)
761 : integer, intent(in) :: nvar ! Number of variables in retVars
762 : character(len=*), intent(in) :: target_name ! Variable being searched for
763 : type(one_dim_read_var), dimension(nvar), intent(in) :: retVars ! Collection being searched
764 :
765 : ! Output Variable
766 : integer :: i
767 :
768 : ! Local Variable(s)
769 : logical :: l_found
770 :
771 : !----------------BEGIN CODE------------------
772 :
773 0 : l_found = .false.
774 :
775 0 : i = 0
776 0 : do while ( i < nvar .and. .not. l_found )
777 0 : i = i+1
778 0 : if( retVars(i)%name == target_name ) then
779 0 : l_found = .true.
780 : end if
781 : end do
782 :
783 0 : if( .not. l_found ) then
784 0 : i = -1
785 : end if
786 :
787 : return
788 :
789 : end function get_target_index
790 :
791 : !=============================================================================
792 0 : function count_columns( iunit, filename ) result( nCols )
793 : ! Description:
794 : ! This function counts the number of columns in a file, assuming that the
795 : ! first line of the file contains only column headers. (Comments are OK)
796 :
797 : ! References:
798 : ! None
799 :
800 : ! Created by Cavyn, July 2010
801 : !-----------------------------------------------------------------------------
802 :
803 : implicit none
804 :
805 : ! External
806 : intrinsic :: index, trim, size
807 :
808 : ! Input Variables
809 : integer, intent(in) :: iunit ! I/O unit
810 : character(len=*), intent(in) :: filename ! Name of the file being read from
811 :
812 : ! Output Variable
813 : integer :: nCols ! The number of data columns in the selected file
814 :
815 : ! Local Variables
816 : integer :: i, k ! Loop Counter
817 : character(len=200) :: tmp ! Temporary char buffer
818 : character(len=200), dimension(50) :: colArray ! Max of 50 columns
819 : logical :: isComment
820 : integer :: status_var ! IO status for read statement
821 :
822 :
823 : ! -------------------------BEGIN CODE-------------------------------------
824 :
825 0 : isComment = .true.
826 :
827 0 : open(unit=iunit, file=trim( filename ), status = 'old' )
828 :
829 : ! Skip all the comments at the top of the file
830 : do while(isComment)
831 0 : read(iunit,fmt='(A)') tmp
832 0 : k = index( tmp, "!" )
833 0 : isComment = .false.
834 0 : if(k > 0) then
835 : isComment = .true.
836 : end if
837 : end do
838 :
839 : ! Go back to the line that wasn't a comment.
840 0 : backspace(iunit)
841 :
842 : ! Count the number of columns
843 0 : nCols = 0
844 0 : colArray = ""
845 0 : read(iunit,fmt='(A)',iostat=status_var) tmp
846 : ! Only continue if there was no IO error or end of data
847 0 : if( status_var == 0 ) then
848 : ! Move all words into an array
849 0 : read(tmp,*,iostat=status_var) (colArray(i), i=1,size( colArray ))
850 :
851 0 : else if ( status_var > 0 ) then
852 : ! Handle the case where we have an error before the EOF marker is found
853 0 : error stop "Fatal error reading data in time_dependent_input function count_columns"
854 :
855 : end if
856 :
857 0 : do i=1,size( colArray )
858 0 : if( colArray(i) /= "" ) then ! Increment number of columns until array is blank
859 0 : nCols = nCols+1
860 : end if
861 : end do
862 :
863 0 : close(iunit)
864 :
865 0 : end function count_columns
866 :
867 : !------------------------------------------------------------------------------
868 0 : end module input_reader
|