Fortran programming help arrays and arithmetic

[twomsta]

Hi folks,

To begin with, I have no past programming experience and have just begun to teach myself programming in FORTRAN 95 and I've hit a wall. I'd be very grateful for any assistance here.

I have around 150 text files with three columns of data (I have attached one as an example, and included as my code so far). Essentially, my aim is to normalize all the precipitation data values in the far right column. But to begin with I have written some code to learn the basics and familiarize myself with programming. So far the code counts how many rows of data are in the text file before using that number to limit how many rows are read (I believe so anyways... as I was having issues with the output when FORTRAN tried to read beyond the length of the text file). No doubt my code so far is unconventional as it was a bit of guessing and modifying, but so far it has done what I wanted.

My question comes down to this:
How exactly do I work with the data values read in as 'prcp', in such a way that I can calculate the mean, standard deviation and eventually the z-score. I assume I need to define an array or dimensions of some sort?

I have read and attempted many tutorials, but has left me a bit lost...Thanks for any help!

program readmet
implicit none 
character(6):: date !I only need the year and month
integer :: endd, period, Nlines
real :: prcp
print*, 'Reads in monthly rainfall values from bureau station data'

!Now calculate the number of rows located in text file (nlines)
Nlines=0
open (5, file='Tprcphq.097091.month.txt')
do
  read (5,*,END=10)
  nlines=nlines+1
  end do
10 close (5)

open (6,file='metout.txt')
write(6,*)' YYYYMM, precipitation (mm)'
write (6,*)'______________________________________________________'

open (5,file='Tprcphq.097091.month.txt')
read (5,*)!read header  
do period= 2, nlines 
  read (5,*)date,endd, prcp
  write (6,*) date, prcp
  print *, date, prcp
     end do
end program readmet

 

[SteamKing]

You don't have to define arrays to calculate the mean.

How would you calculate the mean of a set of numbers if you were doing it by hand?

 

[twomsta]

Thanks for the reply, but unfortunately I don't see what you are getting at ...

Let me put it this way... my output results in two columns of 1277 rows (i.e. the mean would equal the sum off all values in the second column divided by 1277), yet all the values in the second column are read and represented solely by 'prcp'.

Sorry I'm feeling particularly 'slow' today, but I just don't understand how to isolate and identify all the individual values in the second column in such a way that I can carry out simple calculations.

 

[SteamKing]

Yes, I can see that you are 'slow' today.

What is the mean of a string of numbers? It's the sum of all the numbers divided by the number of numbers.

Don't you think you can write some Fortran code which will add up the numbers you read from the data files, count how many numbers you've read, and then calculate what the mean value is? You don't care about storing each individual number for this calculation: you are only interested in the sum and how many numbers there were.

 

[twomsta]

Listen, I have just started teaching Fortran to myself in my spare time, and my mathematical skills are sound outside of programming. I'm not sure if you are trying to come across as demeaning, but I don't need someone to explain the concept of the 'mean' of a set of numbers.

I'm Sure that I am looking at it the wrong way, but if I was able to work with the data that I have read in I wouldn't of bothered posting in this forum.

 

[jbriggs444]

Do you understand how one can update a running total?

 

[DrDu]

Of course you could use arrays, but this isn't necessary here.
E.g. to calculate the mean you could do

prcpmean=0.0e0
do period= 2, nlines
read (5,*)date,endd, prcp
prcpmean=prcpmean+prcp
end do
prcpmean=prcpmean/(nlines-1)
write(*,*) prcpmean

and similarly for the variance etc.

 

[.Scott]

You can also eliminate that first loop.

program readmet
implicit none 
character(6):: date !I only need the year and month
integer :: endd, period, prcpcount
real :: prcp, prcpsum, prcpmean, prcpsquares, prcpvar

print*, 'Reads in monthly rainfall values from bureau station data'

open (6,file='metout.txt')
write(6,*)' YYYYMM, precipitation (mm)'
write (6,*)'______________________________________________________'

open (5,file='Tprcphq.097091.month.txt')
prcpcount = 0
prcpsum=0.0e0
prcpsquares=0.0e0

read (5,*,end=50)!read header
do
  read (5,*,end=50) date, endd, prcp
  write (6,*) date, prcp
  print *, date, prcp
  prcpcount = prcpcount + 1
  prcpsum = prcpsum + prcp
  prcpsquares = prcp*prcp
  end do
50 close (5)

if (prcpcount>1) then
  prcpmean = prcpsum / prcpcount
  prcpvar = (prcpsquares - (prcpsum*prcpsum)/prcpcount)/(prcpcount - 1)
  print *, prcpmean, prcpvar
else
  print *, "Not enough data"
end if

end program readmet

However, now that I've shown you the single pass method of computing the sample variance, I'm going to tell you that your original instincts were good, and (with Fortran) reading through 'Tprcphq.097091.month.txt' twice is the way to go. This is because the 1-pass method of computing the sample deviation can overburden the floating point precision with intermediate results that are quite large values.
Not only that, but you may very well want to hold those values in memory so that you can compute things like the median value.
So...

program readmet
implicit none 
character(6):: date !I only need the year and month
integer :: endd, period, prcpcount, nprcp
real :: prcp, prcpsum, prcpmean, prcpsq, prcpsquares, prcpvar
real, allocatable :: prcps

print*, 'Reads in monthly rainfall values from bureau station data'

open (6,file='metout.txt')
write(6,*)' YYYYMM, precipitation (mm)'
write (6,*)'______________________________________________________'

open (5,file='Tprcphq.097091.month.txt')
prcpcount = 0
prcpsum=0.0e0
prcpsquares=0.0e0

read (5,*,end=50)!read header
do
  read (5,*,end=50) date, endd, prcp
  prcpcount = prcpcount + 1
  prcpsum = prcpsum + prcp
  end do
50 close (5)

if (prcpcount>1) then
  allocate (prcps(prcpcount))
  prcpmean = prcpsum / prcpcount
  rewind 5
  read (5,*)!read header
  do nprcp= 1, prcpcount
    read (5,*) date, endd, prcp
    prcpsq = prcp-prcpmean
    prcpsquares = prcpsq*prcpsq
    prcps(nprcp) = prcp
    write (6,*) date, prcp
    print *, date, prcp
    end do
  prcpvar = prcpsquares/(prcpcount - 1)
  print *, prcpmean, prcpvar
  ... whatever wlse you want to do with the values...
  DEALLOCATE(prcps)
else
  print *, "Not enough data"
end if

end program readmet

One final warning. It has been a couple of decades since I have done any serious Fortran programming, so expect errors in my code.

 

[FactChecker]

Yes, I can see that you are 'slow' today.

This is unnecessarily harsh. Since this is a FORTRAN programming question, let's assume the OP also wants to learn some programming and not just a one-trick solution to a statistics calculation. Ok?

 

[AlephZero]

This is unnecessarily harsh. Since this is a FORTRAN programming question, let's assume the OP also wants to learn some programming and not just a one-trick solution to a statistics calculation. Ok?

In that case, maybe you should give the OP some constructive advice about programming, instead of just criticizing other people's efforts.

 

[FactChecker]

In that case, maybe you should give the OP some constructive advice about programming, instead of just criticizing other people's efforts.

There is no reason to think that the OP ever came back after SteamKing's "slow" comment except to tell off SteamKing. I expect he may not return. There is very little to defend about SteamKing's snide "effort", but if you want to, Ok. Maybe I'm just being slow today.

 

[berkeman]

This is unnecessarily harsh. Since this is a FORTRAN programming question, let's assume the OP also wants to learn some programming and not just a one-trick solution to a statistics calculation. Ok?

You did see that the OP used that term first, right? I don't think SK intended to be insulting. He was just acknowledging what the OP said, before proceeding to provide more help to the OP.

 

[jtbell]

Another point: here on PF, we tend to frown on handing out solutions to problems directly. For the OP's original question of how to calculate the mean without using an array, I would have given a simple example of using a variable to accumulate a sum, and then left it to him to adapt it to his particular situation. Something like this:

sum = 0
n = 0
do
    read (5,*,end=50) data
    sum = sum + data
    n = n + 1
end do
50  print *, n, ' numbers read, with a total of ', sum

 

[twomsta]

Thankyou very much DrDu and Scott, that was what I was after... I had tried writing a similar bit of code earlier on but evidently I was making a mistake somewhere along the line.

To the other comments, there was no need to question my understanding of simple statistical terms in the first place and it wasn't just a one off solution to the mean that I was even after. I wasn't looking for a handout or answer to all my problems, just an example such as jtbell provided. As in the tutorials that I have been doing I have found that examples of similar solutions have been the best way of learning, which I then adapt it to suit my data.

Thanks again for all your help and advice! I should be able to calculate the rest of it.