Strange indexing in Fortran Code

[mushi]

Hi all,

I am totally new to fortran. I have to convert a code that is actually written in the fortran to matlab. There is only one subroutine that is totally mysterious for me. I am posting the part of code I have problem with, below:


1. subroutine cciri(xMODIP,mth,UT,R12,alat,along,foF2,M3000)
2. implicit real*8 (a-h,o-z)
3. real*8 M3000
4. character*10 filena
5. character*2 cm
6. dimension FF0(988),xm0(441),F2(13,76,2),FM3(9,49,2)
7. integer QM(7),QF(9)

8. save

9. data QF/11,11,8,4,1,0,0,0,0/,QM/6,7,5,2,1,0,0/
10. data montha,monthb,Rga/13,14,-10.0D0/

11. if (mth.ne.montha) then
12. write(cm,'(I2.2)') mth+10
13. filena='ccir'//cm//'.asc'
14. open(77,file=filena,status='OLD')

15. read(77,'(4E16.8)') F2,FM3

16. close(77)
17. montha=mth
18. endif

19. if (R12.ne.Rga.or.mth.ne.monthb) then
20. RR2=R12/100.0D0
21. RR1=1.0D0-RR2

22. do i=1,76
23. do j=1,13
24. k=j+13*(i-1)
25. FF0(k)=F2(j,i,1)*RR1+F2(j,i,2)*RR2
26. enddo
27. enddo

Just for reference, here are few lines taken from the files ccirXX.asc

0.56421361E+01 -0.81572428E-01 0.20138618E+00 0.42779669E-01
0.44377025E-01 0.38025160E-02 -0.87150786E-03 0.21708569E-01
-0.23025220E-01 -0.33594195E-02 0.21446483E-03 -0.18263943E-02
-0.12104558E-01 -0.10026588E+01 -0.77602589E+00 0.80115181E+00
-0.27182043E+00 0.54614949E+00 -0.18997940E+00 0.82413152E-01
-0.38998052E-01 -0.97885408E-01 0.13837354E+00 0.10652158E+00

Can anybody please explain lines 15 & 25 to me. I mean on line 15, variables F2 and FM3 are loaded as variables while on line 25 they are treated as 3D arrays? How is this possible? What is the logic behind this code? Mind you that the code is running absolutely fine in frotran.
Please help!.

 

[Well_InTheory]

6. dimension FF0(988),xm0(441),F2(13,76,2),FM3(9,49,2)

declares the variables F2 and FM3 to be arrays with the listed dimensions.

15. read(77,'(4E16.8)') F2,FM3

reads the values into the variables according to the style

(4E16.8)

i.e. 4 blocks of E16.8 data, as you see in the data file. The read statement is ambiguous (I'd have used an iterator in a loop) but I'd guess that it will fill the three F2 slots then the first FM3 slot... ?

My '77 is a little rusty, but without a do-loop around that read, or iterators, I'm not certain that the data will be filled correctly. You could try printing to screen F2 and FM3 directly after the read and seeing for yourself what is happening. I'd be worried about the '2' in the dimensions of each - are they being filled sequentially?

 

[mushi]

Here is the remaining part of the subroutine:

28. do i=1,49
29. do j=1,9
30. k=j+9*(i-1)
31. xm0(k)=FM3(j,i,1)*RR1+FM3(j,i,2)*RR2
32. enddo
33. enddo

34. Rga=R12
35. monthb=mth
36. endif

I am showing it because this makes the situation more complex in a sense that line 31 is similar to line 25 in the previous part but this time it is FM3.

@Well_InTheory: Yes bro that '2' is the main problem. I can't understand why the author put it there. That '2' is actually making it 3D matrix.

 

[Well_InTheory]

The fact that it's 3D isn't a problem - that last dimension just corresponds to whether the variable will multiply RR1 or RR2, my problem is that given the do-loops, ALL of the elements should be populated, but there's complete ambiguity over which will get filled when.

Note: as long as the read statement is filling these variables the way the author intended, there's nothing wrong here. The dimension 3 arrays could also have been made as dimension 2 arrays;

FM3_1(1:9,1:49) which multiplies RR1

and

FM3_2(1:9,1:49) which multiplies RR2

which would also have been easier to read in with less ambiguity. In any case, this is what you have, so print the values after the read and map out the data file yourself, then you'll know what goes where.

 

[mushi]

Thanks Well_InTheory for your time buddy.

Please explain what do u mean by:

'The fact that it's 3D isn't a problem - that last dimension just corresponds to whether the variable will multiply RR1 or RR2'

I am trying to convert this script in MATLAB and that's where they appear to be 3D. Kindly explain your statement above.

 

[jim hardy]

can't be this simple, can it?


Had you not said you are "totally new to Fortran" i'd remain silent because I'm totally rusty on it.

I think a Fortran "read" of an already dimensioned variable will try to fill the whole array.

so line 6 tells line 15 how much reading he has to do.

my last Fortran was about 1986 - maybe somebody who's fresher can confirm or correct.
i would have used at least an implied do , for sake of next guy, or commented.

old jim

 

[mushi]

Thanks Jim for the reply.

This subroutine has not been written by me. I got it from my tutor to convert it and use it for my own assignment in matlab. That's the reason why I am totally confused by the dimension set for the F2 and FM3 arrays and then their strange 3D use on lines 25 and 31. Well_InTheory has given me some good points to work on. But still I am not clear on his last post and waiting for his response.

Thanks anyways.
Keep posting mate.

 

[jim hardy]

i hope I'm on the beam here, if not just dismiss it.

here's a link that seems appropos ( but i may not even understand the question. )

http://www.idris.fr/data/cours/lang/fortran/f90/F77.html#p7.4

see section 7.5, multidimensional arrays

The unsubscripted name of a multi-dimensional array can be used, like that of a one-dimensional array, in input/output and DATA statements to refer to all its elements, but it is essential to know their order. The elements are referenced in the order of their positions in the computer's memory. For a one-dimensional array, the elements occur, as we might expect, in increasing order of their subscripts, but for multi-dimensional arrays, the ordering is less obvious. The rule is that the elements are ordered with the first subscript increasing most rapidly, then the next and so on, the last subscript increasing most slowly. Thus if MARK is declared as:
INTEGER MARK(5,100)

its elements are ordered in memory as shown above, and the statement:
READ *,MARK

is equivalent to:
READ *, ((MARK(I,J),I=1,5),J=1,100)

hmmm - indeed 25 and 31 make 3-d's into one-d's.

Goood Luck !

And, thanks for your kind words.

 

[mushi]

Its not about the normal two dimensions, its about how the variables/arrays F2 and FM3 are loaded from the file on line 15 and read on lines 25 and 31 in this particular code. For me, they were loaded as simply 2D arrays but read as 3D arrays. I am so confused because I am used to work in Matlab where you have to be consistent with the dimensions of the matrix throughout the code.

Can you please simply tell me how can I reproduce the loading (line 15) and and reading (lines 25,31) in matlab.

I hope I a clear this time.

 

[AlephZero]

The read statement is ambiguous (I'd have used an iterator in a loop) but I'd guess that it will fill the three F2 slots then the first FM3 slot... ?

My '77 is a little rusty, but without a do-loop around that read, or iterators, I'm not certain that the data will be filled correctly.

There is nothing "ambiguous" about the read statement and that explanation is wrong.

It reads values into every element of F2, then every element of F3, taking 4 values from each line of the file and reading as many lines from the file as required.

The array elements are read in order
F2(1,1,1) F2(2,1,1) ... F2(13,1,1)
F2(1,2,1) F2(2,2,1) ... F2(13,76,1)
F2(1,1,2) F2(2,1.2) ... F2(13,76,2)
and then the elements of F3.

 

[gsal]

According to the instructions quoted by Jim, sounds like you could replace the code, more explicitely with

do i = 1,13
do j = 1,76
do k = 1,2
   read(77,*) F2(i,j,k)
end do
end do
end do

do i = 1,9
do j = 1,49
do k = 1,2
   read(77,*) FM3(i,j,k)
end do
end do
end do

If you have the ability to compile fortran, you could write a 20-line program to test this. If you do not have the ability to compile fortran...acquire it! There is a nice free compiler called g95, you could download for Linux or Windows.

 

[mushi]

Thanks AlephZero and gsal. I think gsal made my life easy by putting some coding thing which I think is the best explanation of what others are talking about. I will definitely try gsal way. I hope this will what I thought it would be after reading others' comments.

keep posting guys.
Thanks

 

[AlephZero]

According to the instructions quoted by Jim, sounds like you could
replace the code, more explicitely with

do i = 1,13
do j = 1,76
do k = 1,2
   read(77,*) F2(i,j,k)
end do
end do
end do
(etc)

Sorry, but that is wrong. It should be

do k = 1,2
do j = 1,76
do i = 1,13
   read(77,*) F2(i,j,k)
end do
end do
end do
(etc)

But that version still won't work properly unless you change the data file. The original code reads 4 numbers from each line of the file, but this code will only read the first 2 numbers from each line, because each READ statement starts reading the next line of the file. It will eventually crash, when it tries to read past the end of the file.

There is nothing wrong with the original code, and you don't need to change it. The only "problem" is that you didn't understand what it does.

 

[mushi]

Dear AlephZero, it is not the matter of just understanding the code, if you read the full post, it is actually to convert the code into MATLAB so that i may be able to use it for my assignment which is actually in matlab. I will really appreciate if you could suggest something in this regard.
Thanks anyways.

 

[gsal]

o.k., I admit it...it was 2 minutes before a meeting when I posted the nested loops...not that I would have caught the fact that I had my order reverse...good catch.

I think that if the read statement gets enhanced with a format string it would work just fine; it may be right that as it is it will not work ...with the list directed input ( the asterisk ' * ') , but if you write something like

read(77, '(E16.8)')

then it will work...I think fortran will simply just read one number at a time without going to the next line as indicate before...

...like I said...if you could acquire the capability to compile fortran (easy!) you could be verifying what the code does here and there and carry out a better translation/port of the code.

...by the way, did you know that Matlab is written in Fortran? I don't use Matlab, so I don't know; but I wonder if you have the ability to use user-typed-compiled fortran programs as plug-in functions and so you wouldn't have to be porting anything and including the possibility to make a mistake due to lack of understanding of the original code...just a thought...I probably wouldn't recommend doing input/output at the fortran level, though, just number crunching.

 

[Well_InTheory]

There is nothing "ambiguous" about the read statement and that explanation is wrong.

It reads values into every element of F2, then every element of F3, taking 4 values from each line of the file and reading as many lines from the file as required.

The array elements are read in order
F2(1,1,1) F2(2,1,1) ... F2(13,1,1)
F2(1,2,1) F2(2,2,1) ... F2(13,76,1)
F2(1,1,2) F2(2,1.2) ... F2(13,76,2)
and then the elements of F3.

Well, it *IS* ambiguous to read if you don't inherently know the rule it will go by, made worse that the choice to have 4 datum per line is thus entirely arbitrary. Especially if you're more used to an iterated approach which might fill a series of 1D arrays where the data progresses down a file. As presented it would be rather impossible to look at the data file and see where the dimension changes - it may as well have been written as a 1D array.

My explanation was a guess, and was cautioned as such.

Its not about the normal two dimensions, its about how the variables/arrays F2 and FM3 are loaded from the file on line 15 and read on lines 25 and 31 in this particular code. For me, they were loaded as simply 2D arrays but read as 3D arrays. I am so confused because I am used to work in Matlab where you have to be consistent with the dimensions of the matrix throughout the code.

Can you please simply tell me how can I reproduce the loading (line 15) and and reading (lines 25,31) in matlab.

I hope I a clear this time.

I think this is where the trouble arose - why do you think they are loaded as 2D? They are declared as 3D in the dimension statement. This is entirely consistent.

 

[mushi]

In my view 2D corresponds to (row, column). I found the data in the file in the same manner i.e. in rows and columns. How would you define a 3rd dimension in that data file no matter how the arrays were defined in the dimension statement?

PS: This is not my question this is actually matlab's question and I am trying to answer matlab.

 

[barrowes]

Here is a little MATLAB code and an accompanying fortran snippet to help with the translation:

First the MATLAB code:
F2=[1:13*76*2];

fid=fopen('test1','w');
fprintf(fid,'%16.8f %16.8f %16.8f %16.8f \n',F2);
fclose(fid)

f2=reshape(F2,[13,76,2]);
f2(1:13,1,1)
f2(1,1:10,1)
f2(1,1,1:2)


fid=fopen('test1','r');
F2r=fscanf(fid,'%f');
f2r=reshape(F2r,[13,76,2]);
f2r(1:13,1,1)
f2r(1,1:10,1)
f2r(1,1,1:2)

fclose(fid)





And the fortran:

program test1f
real F2(13,76,2)
open(77,file='test1',status='OLD')
read(77,'(4E16.8)') F2
print *,'F2(1:13,1,1)=',F2(1:13,1,1)
print *,'F2(1,1:10,1)=',F2(1,1:10,1)
print *,'F2(1,1,1:2)=',F2(1,1,1:2)
close(77)


end program test1f



HTH. Also you could try f2matlab for future conversion projects.

Ben Barrowes

 

[gsal]

mushi:

I lost you...I don't even know why we are talking 2D...maybe somebody earlier said that the matrices are being loaded in 2D...but that is not true...the matrices were declared as 3D matrices and are being loaded one item at a time, sequentially in "Fortran" order, the numbers themselves could have been 1 per row, 2, 3, or 4 per row, or ALL in a single row.

Needless to say, a file looks to us like a piece of paper in which you can write numbers in rows...if you happen to have a 2D matrix, well, you could be friendly enough to arrange the numbers on the file so that they actually reflect your matrix...but if you just dump it from Fortran, it would look one way and if you dump it from C, it would look some other way.

What if you have a 3D matrix? or a 4D matrix or a 5D matrix? How do you write it to a file? Well...take your pick! You could do anything you understand.

You could go the human route and explicitly write one entry per row clearly specifying the index the entry goes into. Or so many other ways..it would be a matter of taste.

Or, you could go the computer route and simply 'dump' the matrix into a file and then let it be read in the same way. Fortran has this ability of just writing an entire matrix into a file and reading it back in. When I studied java, I learned about serialization where you could write to file any object that you had, could be a matrix...then again, in order to read it back you need to know its class in advance. In python with 'pickle', you can accomplish the same thing, but nicely enough, I don't think you need to know in advance what you are reading.

Almost lastly, I am not quite sure what you are asking in your last question...can you be more clear?

Lastly...can you declare a 3D matrix in MATLAB and write it out and see what MATLAB does? This could be helpful.

 

[jim hardy]

Please indulge me a simplification here.

I think we over-anthropomorphise our computers.

Put yourself in the shoes of the computer's ALU.

There is no such thing as dimensions, just memory addresses and those are sequential.
Dimensions exist only in the mind of the programmer.

The ALU loads numbers from the input file, in sequence, and stores them in physical memory , i would guess also in sequence.. if it's DMA they're certainly both in sequence, on any machine i ever saw.

The programmer has chosen to refer in his mind to the location where those numbers are stored not by the sequence they were encountered in the input file, but as if they were arranged in a nice orderly multidimensional matrix.
It is the job of the compiler(? assembler?) to translate between the physical memory location needed by ALU and the programmer's imagined location.

In your case* he imagined 13 rows, 76 columns and 2 layers. That's 1976 entries, interestingly twice the size of array FFO where he writes results later.
(* pardon me for not knowing exactly what Fortran does - i might have order of rows columns and layers not proper, )

let me oversimplify for illustration - take a 10X10X10 matrix.
the first ten numbers from a file might fill one row of top layer,
and the next 90 numbers fill the remaining nine rows of that layer,
and the next 900 numbers fill the remaining nine layers,
in the mind of the programmer.
But they are stored unidimensionally in memory.

Again, it is the job of the compiler to translate between programmer's imagined location and the physical memory location to be asserted by ALU when it loads..

To get at the next to last entry the ALU must retrieve the 999th entry,
and whether it got there by a sequential count,
or by a calculation something like: layer*100 + row*10 + column
is not visible to us as Fortran programmers.


13 rows, 76 columns and 2 layers... each layer is 13 X 76 = 988 entries;;
interestingly do i read that lines 22-27 make a one dimensional list out of the products of layer 1 and layer 2, vertically aligned entry by entry, resulting in a 988 entry one dimension list?
Could he have as easily read the file into a one dimensioned 1976 entry array, called it Fred,
then
For i=1 to 988
product(i) = Fred(i) * Fred(i + 988)
next i
?


Please excuse these ramblings of an old guy. Its just i sense confusion, and the remedy for that is always back to basics.



Maybe this will help a beginner...
If i am pestering you just say so.



old jim

 

[Well_InTheory]

In my view 2D corresponds to (row, column). I found the data in the file in the same manner i.e. in rows and columns. How would you define a 3rd dimension in that data file no matter how the arrays were defined in the dimension statement?

PS: This is not my question this is actually matlab's question and I am trying to answer matlab.

The formatting of the data file is irrelevant - the entries could all be in a single column for all Fortran cares, as long as it knows how to read it.

Since the third index only runs over [1,2] you could simply have two blocks of code for which in each you hold that last index fixed.

The simpler option would be to not use a 3D array where you simply want to use the two indices separately, and have

F2_1 and F2_2, FM3_1 and FM3_2 where the _n corresponds to which will multiply RRn, then you can have your (row,column) data file map nicely onto the actual entries.

Of course, this mapping will not be what is in the data file now, you'll need to reformat it. The easiest way would be to have the current program output the arrays as you will want them in Matlab and create a new datafile.