Interpolating & Integrating Data from File to 2D Array

[VoidSerpent]

I have this code that i have been working on that will interpolate and integrate some data. I need to read this data from a file with unknown rows and columns to a two dimensional array/matrix so that the data can be broken down further into one dimensional arrays from the parent 2D array's columns. so a quick example is the data may look like this in the file:
1 2 3
4 5 6
7 8 9
and so i would read it into an allocatable matrix becoming [3,3]. That is where my programming logic breaks down with FORTRAN. I have found codes and literature on allocatable 1D arrays and some on reading files to 2D arrays but not allocatable 2d. so my thought is to obtain an array/matrix that at point m(2,2) = 5 i need a code roughly like this:

	   ! reading file to an array/matrix
      do h=1, number_rows
         do k=1, number_columns
           read(1, *, iostat = io) a
         enddo
      enddo

but again number_rows and number_columns are unknown and if i included an if or do loop that ended at the end of the line of each row i don't know what might be needed for that sort of declaration. because it would be nice to have a nested do loop that the inside one increments all the way across a row till the end of the line then the outer one will increment so that the program goes to a new row/line so that the inner loop increments again repeating till end of file. were it a file that i could put into a 1D array the code would be like so:

DO
READ(1, *, IOSTAT = IO) A
IF (IO < 0 ) EXIT
N = N + 1
ALLOCATE( OldX( SIZE(X) ) )
OldX = X ! entire array can be assigned
DEALLOCATE( X )
ALLOCATE( X(N) )
X = OldX
X(N) = A
DEALLOCATE( OldX )
END DO

so my question to anyone that knows FORTRAN and understands what i am trying to accomplish is how do i read a 2D file into an allocatable 2D array? If it helps clarify anything i am also attaching the rest of the code before and after where this part is needed.

      Program main
      implicit none      
      real, dimension(:,:), allocatable :: mc, mr, oldm
      real a
      integer io, nc, nr, h, k
      character(30) :: filename
	  real, dimension(:,:), allocatable :: alt, temp, nue, oxy
      integer locationa, locationt, locationn, locationo, i
      integer nend
	  real dz, z, integral
	  real alti, tempi, nuei, oxyi
	  integer y, j

      allocate( m(0, 0) ) ! size zero to start with?

      nn = 0
      j = 0

       write(*,*) 'Enter input file name: '

       read(*,*) filename

       open( 1, file = filename )
	   !this is the needed code
	   !
	   ! reading file to an array/matrix
      do h=1, nr
         do k=1, nc
           read(1, *, iostat = io) a
         enddo
      enddo
	  !
	  !
	  !between these comments
	  ! Decompose matrix array m into column arrays [1,n]

      write(*,*) 'Enter Column Number for Altitude'
	  read(*,*) locationa
	  write(*,*) 'Enter Column Number for Temperature'
	  read(*,*) locationt
	  write(*,*) 'Enter Column Number for Nuetral Density'
	  read(*,*) locationn 
	  write(*,*) 'Enter Column Number for Oxygen density'
      read(*,*) locationo

      nend = size(m, locationa) !length of column #locationa, or any column from the data file
	  
      do i = 1, nend

       alt(i, 1) = m(i, locationa)

       temp(i, 1) = log(m(i, locationt))

       nue(i, 1) = log(m(i, locationn))

       oxy(i, 1) = log(m(i, locationo))

      enddo



      ! Interpolate Column arrays, Constant X value will be array ALT with the 3 other arrays

      !real dz = size(alt)/100, z, integral = 0
	  !real alti, tempi, nuei, oxyi
	  !integer y, j = 0
	  dz = size(alt)/100


      do z = 1, 100, dz
	  y = z !with chopped rounding alt(y) will always be lowest integer for smooth transition.
	  alti = alt(y, 1) + j*dz ! the addition of j*dz's allow for all values not in the array between two points of the array. 

       tempi = exp(linear_interpolation(alt, temp, size(alt), alti))

       nuei = exp(linear_interpolation(alt, nue, size(alt), alti))

       oxyi = exp(linear_interpolation(alt, oxy, size(alt), alti))
	   j = j + 1



       !Integration

       integral = integral + tempi*nuei*oxyi*dz 

      enddo


      end program main


      !Functions

      real function linear_interpolation(x, y, n, x0)

       implicit none

       integer :: n, i, k

       real :: x(n), y(n), x0, y0

       k = 0


      do i = 1, n-1

       if ((x0 >= x(i)) .and. (x0 <= x(i+1))) then  

        k = i ! k is the index where: x(k) <= x <= x(k+1)
        exit ! exit loop

       end if

      enddo


      if (k > 0) then  ! compute the interpolated value for a point not in the array

       y0 = y(k) + (y(k+1)-y(k))/(x(k+1)-x(k))*(x0-x(k))

      else

       write(*,*)'Error computing the interpolation !'

       write(*,*) 'x0 =',x0, ' is out of range <', x(1),',',x(n),'>'

      end if

      ! return value

         linear_interpolation = y0

      end function linear_interpolation

 

[rcgldr]

I'm not that familiar with current versions of Fortran. Is is possible to obtain the size of the file, then to allocate a buffer and read all of the file into the buffer as Ascii data, then scan the buffer to determine the number of rows and columns, then allocate the matrx and scan the buffer again, converting the strings into numbers in the matrix? Or after scanning the buffer to determine matrix size, allocate the matrix, then close and reopen the file, but on the second pass, read each line as a row of numbers to input into the rows of the matrix?

 

[AlephZero]

The easy way to do this would be to change the file, so the first line contains the number of rows and columns.

Then you can read the size of the array from the first line, allocate the array, then read the all the data with one READ statement in free format. That will read as many lines from the file as it needs, to get all the numbers.

If you can't change the file format, and you are assuming each line of the file corresponds to one row of the matrix, I would do it by reading the file twice.

First read each line as a character string, count the number of data items (i.e. strings separated by blanks) on the lines, and count the number of lines. When you get to the end of file, Alllocate space for the data, rewind the file, and read the numbers into the array you allocated.

 

[VoidSerpent]

@rcgldr, i researched inquiry and length/size code tid bits and none exist to determine the file until it is read. I found some things that suggest 'rewind' as a possible tool to use but nothing on applying it toward columns, closest i found was for the number of rows with

DO
READ(13,*, iostat = error)
IF (error == -1) EXIT
n = n + 1
END DO

which maybe i will utilize for the rows which are far more numerous than the columns anyway in the data i am currently working with. But if i start adding variable columns it would be nice to not have it user defined.

@AlphaZero do you have a link to how i would read the file into strings? I will scan my books and pdfs now but if you have a good reference site i would appreciate it, obviously some of my literature doesn't cover the things i need.

 

[rcgldr]

@rcgldr, i researched inquiry and length/size code tid bits and none exist to determine the file until it is read.

One option if accessable in Fortran, is to open a file for read, seek to end of file, then get the file position, which would be the number of bytes (or number of bytes -1), then do a rewind (seek to beginning of file), then scan the file for the number of rows. You didn't mention what OS (Linux, Windows, ... ) you're doing this on, or if you're trying to keep the code OS independent.

 

[VoidSerpent]

Thanks for the assist on reading in the data, i believe i have a code that will do it properly now but i am still getting compilation errors. all of them refer to my subroutine they are

real subroutine linear_interpolation(x, y, n, x0)
1
Error: Unexpected data declaration statement at (1)
int.f:125.72:

integer :: n, b, k
1
Error: Unexpected data declaration statement at (1)
int.f:127.72:

real :: x(n,1), y(n,1), x0, y0
1
Error: Unexpected data declaration statement at (1)
int.f:158.29:

linear_interpolation = y0
1
Error: 'linear_interpolation' at (1) is not a variable
int.f:160.9:

end subroutine linear_interpolation
1
Error: Expecting END PROGRAM statement at (1)
int.f:132.9:

do b = 1, n-1
1
Warning: Deleted feature: Loop variable at (1) must be integer

at first i thought the unexpected declarations may have been the arrays passing into my subroutine were not equitable to the arrays used in it, but it appears i fixed all that and am still getting those errors. As for the "linear_interpolation" not being a variable, i don't know how to declare it any other way to be able to return my desired value. and i really have no clue what it is asking for with the "end program" statement. my new code is as follows

 Program main
	 ! implicit none      
     
      real, dimension(:,:), allocatable :: m
      character(100) columncount
      integer io, nc, nr
      character(30) :: filename
      real, dimension(:,:), allocatable :: alt, temp, nue, oxy
      integer locationa, locationt, locationn, locationo, i
      integer h, g, error, int
      real dz, z, integral
      real alti, tempi, nuei, oxyi, t, nu, o
      integer l, j       write(*,*) 'Enter input file name: '
       read(*,*) filename
       open(unit=10,file=filename,access='sequential',form='formatted')

 

      !reading in data file
	  
	  read(1, *) columncount !reads first line to be able to determine # of columns
      rewind(1)
	  
	  
	  nc = len(columncount)
      do while (columncount(nc:nc) .eq. ' ')
        nc = nc - 1
      enddo
	  
	   
	  nr = 0
      do
      read(1,*, iostat = error)
      if (error == -1) exit
      nr = nr + 1
      enddo
      rewind(1)
	  
	  do h=1, ubound (m, 1)
      read (10, '(5I4)')  (m (h, g), g=1, ubound (m, 2))
      end do
      ! Decompose matrix array m into column arrays [1,n]

      write(*,*) 'Enter Column Number for Altitude'
	  read(*,*) locationa
	  write(*,*) 'Enter Column Number for Temperature'
	  read(*,*) locationt
	  write(*,*) 'Enter Column Number for Nuetral Density'
	  read(*,*) locationn 
	  write(*,*) 'Enter Column Number for Oxygen density'
      read(*,*) locationo

      !nr length of column or # of rows
	  
      do i = 1, nr

       alt(i, 1) = m(i, locationa)

       temp(i, 1) = log(m(i, locationt))

       nue(i, 1) = log(m(i, locationn))

       oxy(i, 1) = log(m(i, locationo))

      enddo
      ! Interpolate Column arrays, Constant X value will be array ALT with the 3 other arrays	  dz = size(alt)/100
      z=1

      do int = nint(z), 100
	  l = z !with chopped rounding alt(l) will always be lowest integer for smooth transition.
	  alti = alt(l, 1) + j*dz ! the addition of j*dz's allow for all values not in the array between two points of the array. 
      t = linear_interpolation(alt, temp, nr, alti)
       tempi = exp(t)
      nu = linear_interpolation(alt, nue, nr, alti)
       nuei = exp(nu)
      o = linear_interpolation(alt, oxy, nr, alti)
       oxyi = exp(o)
	   j = j + 1
	   z = z + dz       !Integration

       integral = integral + tempi*nuei*oxyi*dz 

      enddo      !end program main

      !contains
      !Functions

      real subroutine linear_interpolation(x, y, n, x0)

       !implicit none

       integer :: n, b, k

       real :: x(n,1), y(n,1), x0, y0

       k = 0      do b = 1, n-1

       if ((x0 >= x(b,1)) .and. (x0 <= x(b+1,1))) then  

        k = i ! k is the index where: x(k) <= x <= x(k+1)
        exit ! exit loop

       end if

      enddo      if (k > 0) then  ! compute the interpolated value for a point not in the array

       y0 = y(k,1) + (y(k+1,1)-y(k,1))/(x(k+1,1)-x(k,1))*(x0-x(k,1))

      else

       write(*,*)'Error computing the interpolation !'

       write(*,*) 'x0 =',x0, ' is out of range <', x(1,1),',',x(n,1),'>'

      end if

      ! return value

         linear_interpolation = y0
      return
      end subroutine linear_interpolation
      end program main

I didn't know if i should start a new thread for this new problem or not. Any help is appreciated, even if it is just in reference to a better literature source, that is how i have tweaked most of this code already.

EDIT: on another thread Mark44 mentioned that subroutines shouls be outside the main program, mine originally were but then i was ONLY getting
Program main
1
int.f:109.72:

real subroutine linear_interpolation(x, y, n, x0)
2
Error: Two main PROGRAMs at (1) and (2)
so i hope that helps in the diagnosis of my problem

 

[Mark44]

You have declared linear_interpolation() as a subroutine, but you are using it as if it had been declared a function.

[color="red"]t = linear_interpolation(alt, temp, nr, alti)[/color]
tempi = exp(t)
[color="red"]nu = linear_interpolation(alt, nue, nr, alti)[/color]
nuei = exp(nu)
[color="red"]o = linear_interpolation(alt, oxy, nr, alti)[/color]

A subroutine causes something to happen, but does not return a value, so should not be on the right side of an assignment statement or otherwise used as a value.

On the other hand, your subroutine attempts to return a value, so it should be written as a function, not a subroutine.

 

[VoidSerpent]

ok, ill change it back to a function. I had people tell me their codes never use "functions" and to change it to a subroutine so i did. but i got my code complied finally. i need to work on reading the file in and counting the items in the string. Is there a way to only read one line from the file into the string, the numbers are of form '1.123E-12'? i thought this would work but in my test file i have about 13 columns and it returns nc as 4. and when i write the string to the screen its one number that isn't even in the file.

write(*,*) 'Enter input file name: '
read(*,*) filename
open(unit=10,file=filename,access='sequential',form='formatted')

 
!reading in data file
	  
read(10, *) columncount !reads first line to be able to determine # of columns
rewind(10)
	  
	  
nc = len(columncount)
do while (columncount(nc:nc) .eq. '  ')
        nc = nc - 1
enddo

 

[Mark44]

What are you doing here? Can you explain it to me in words?

Also, what you have in the single quotes looks to be two spaces.

nc = len(columncount)
do while (columncount(nc:nc) .eq. '  ')
        nc = nc - 1
enddo

I had people tell me their codes never use "functions" and to change it to a subroutine so i did.

I wonder if it was because they didn't understand the difference between functions and subroutines.

 

[VoidSerpent]

What are you doing here? Can you explain it to me in words?

Also, what you have in the single quotes looks to be two spaces.

nc = len(columncount)
do while (columncount(nc:nc) .eq. '  ')
        nc = nc - 1
enddo

Well i took that from a PDF on Intro to FORTRAN. This is the section;

The blank padding at the end of the string is counted when you use the
LEN() function to find the string's length, or when you WRITE the string.
To find the "true" length of the string use:

integer function strlen(st)
integer i
character st*(*)
i = len(st)
do while (st(i:i) .eq. ' ')
i = i - 1
enddo
strlen = i
return
end

Strings don't come initialized with blanks, if the compiler initializes
them (VMS, Sun) they are initialized to NULs. Note that some terminals
(e.g. VTnnn) ignore NUL characters and if such a string is written to the
screen there will be no visible output (except the start of a new line).

and so i was trying to use that to obtain the number of columns that i will need to properly read in the file to a matrix array. Here is my some what working code, it compiles but now i am running into these run time errors. Thanks in advance and up to now for the help.

      !Program main
	 ! implicit none      
     
      real, dimension(:,:), allocatable :: m
      character(10000) columncount
      integer io, nc, nr
      character(30) :: filename
      real, dimension(:,:), allocatable :: alt, temp, nue, oxy
      integer locationa, locationt, locationn, locationo, i
      integer h, g, error, int
      real dz, z, integral
      real alti, tempi, nuei, oxyi, t, nu, o
      integer l, j


       write(*,*) 'Enter input file name: '
       read(*,*) filename
       open(unit=10,file=filename,access='sequential',form='formatted')

 

      !reading in data file
	  
	  read(10, *) columncount !reads first line to be able to determine # of columns
      rewind(10)
	  
	  
	  nc = len(columncount)
      do while (columncount(nc:nc) .eq. '  ')
        nc = nc - 1
      enddo
	  
	   
	  nr = 0
      do
      read(10,*, iostat = error)
      if (error == -1) exit
      nr = nr + 1
      enddo
      rewind(10)
	  
	  do h=1, ubound (m, 1)
      read (10, '(5I4)')  (m (h, g), g=1, ubound (m, 2))
      end do
      
      write(*,*) columncout
      write(*,*) nc
      write(*,*) nr


      ! Decompose matrix array m into column arrays [1,n]

      write(*,*) 'Enter Column Number for Altitude'
	  read(*,*) locationa
	  write(*,*) 'Enter Column Number for Temperature'
	  read(*,*) locationt
	  write(*,*) 'Enter Column Number for Nuetral Density'
	  read(*,*) locationn 
	  write(*,*) 'Enter Column Number for Oxygen density'
      read(*,*) locationo

      !nr length of column or # of rows
	  
      do i = 1, nr

       alt(i, 1) = m(i, locationa)

       temp(i, 1) = log(m(i, locationt))

       nue(i, 1) = log(m(i, locationn))

       oxy(i, 1) = log(m(i, locationo))

      enddo



      ! Interpolate Column arrays, Constant X value will be array ALT with the 3 other arrays


	  dz = size(alt)/100
      z=1

      do int = nint(z), 100
	  l = z !with chopped rounding alt(l) will always be lowest integer for smooth transition.
	  alti = alt(l, 1) + j*dz ! the addition of j*dz's allow for all values not in the array between two points of the array. 
      t= linear_interpolation(alt, temp, nr, alti)
       tempi = exp(t)
      nu= linear_interpolation(alt, nue, nr, alti)
       nuei = exp(nu)
      o= linear_interpolation(alt, oxy, nr, alti)
       oxyi = exp(o)
	   j = j + 1
	   z = z + dz


       !Integration

       integral = integral + tempi*nuei*oxyi*dz 

      enddo


      end !program main

      !contains
      !Functions




      function linear_interpolation(x, y, n, x0)

       !implicit none

       integer n, b, k

       real x(n,1), y(n,1), x0, y0

       k = 0
       

      do b = 1, n-1

       if ((x0 >= x(b,1)) .and. (x0 <= x(b+1,1))) then  

        k = i ! k is the index where: x(k) <= x <= x(k+1)
        exit ! exit loop

       end if

      enddo


      if (k > 0) then  ! compute the interpolated value for a point not in the array

       y0 = y(k,1) + (y(k+1,1)-y(k,1))/(x(k+1,1)-x(k,1))*(x0-x(k,1))

      else

       write(*,*)'Error computing the interpolation !'

       write(*,*) 'x0 =',x0, ' is out of range <', x(1,1),',',x(n,1),'>'

      end if

      ! return value
      linear_interpolation = y0
      !return
      end !subroutine linear_interpolation

and to reiterate the data from the file is roughly [168,13] and of the form '1.123E-12', i saw that i may need some format declaration for scientific notation, when i have time ill add that to the code where i believe it needs to go.

 

[VoidSerpent]

So i have been working on my code, and i realized that the program is not reading an entire line, just the first letter of each line. I created a bunch of dummy files to test how it was reading. i also got this code

      PROGRAM readfile
      IMPLICIT NONE
      REAL, DIMENSION(:), ALLOCATABLE :: mydatar
      INTEGER, PARAMETER :: maxrecs = 10000
      INTEGER :: J, NR, ios
      CHARACTER(100) :: filename
      CHARACTER(1) :: junk
      write(*,*) 'Enter name of file to read in...'
      read(*,*) filename
      !Determine total number of lines in file
      NR = 0
      OPEN(UNIT=1,FILE=filename)
      DO J=1,maxrecs
      READ(1,*,IOSTAT=ios) junk
      IF (ios /= 0) EXIT
      IF (J == maxrecs) THEN
      write(*,*) 'Error: Maximum number of records exceeded...'
      write(*,*) 'Exiting program now...'
      STOP
      ENDIF
      NR = NR + 1
      ENDDO
      REWIND(1)
      !Now we can allocate data variables
      ALLOCATE(mydatar(NR))
      !Now read data into mydata
      DO J=1,NR
      READ(1,*) mydatar(J)
      write(*,*) mydatar(J)
      ENDDO
      CLOSE(1)
      END PROGRAM readfile

to do a secondary comparison on reading in the file. both lead me to think there is something in the formatting of the files because they both only read the first number in each line. I have saved the files as; .txt, .csv, and .ods as well as changing the field and text delimiters to no avail. Any suggestions?