ccccccc Programa per a buscar la diferencia entre dues matrius
	real*8 matriu1,matriu2,resta
	real*8 resta_min,resta_max,migresta,migquadresta,varresta
	real*8 mig1,migquad1,mig2,migquad2,mig12,migquad12
	real*8 var1,var2,var12,corr,mae,rmse
	integer i,j
	integer n_a1,n_a2,n_as,n_ubic1,n_ubic2,n_ubics
	integer nx,ny,dim1,dim2
	character*80 arxiu_1,arxiu_2,arxiu_sortida
	character*200 ubicacio1,ubicacio2,ubicacio_sortida
	dimension matriu1(300,300),matriu2(300,300),resta(300,300)
	dimension dim1(2),dim2(2)
	
	call llegir(ubicacio1,ubicacio2,arxiu_1,arxiu_2,
     &ubicacio_sortida,arxiu_sortida)
	
	call llargada_String(arxiu_1,80,n_a1)
	call llargada_String(arxiu_2,80,n_a2)
	call llargada_String(arxiu_sortida,80,n_as)
	call llargada_String(ubicacio1,200,n_ubic1)
	call llargada_String(ubicacio2,200,n_ubic2)
	call llargada_String(ubicacio_sortida,200,n_ubics)

	write(6,*)'Tenim les dimensions de les matrius...'
	call dimensions(arxiu_1,ubicacio1,dim1(1),dim1(2))
	write(6,*)'Matriu 1 a arxiu... ',arxiu_1(1:n_a1),dim1(1),
     &' x ',dim1(2)
	write(6,*)' '
	call dimensions(arxiu_2,ubicacio2,dim2(1),dim2(2))
	write(6,*)'Matriu 2 a arxiu... ',arxiu_2(1:n_a2),dim2(1),
     &' x ',dim2(2)
	write(6,*)' '
	
	if (dim1(1).ne.dim2(1)) then
	write(6,*)'Dimensions de les matrius diferents !!!!'
	stop
	endif
	if (dim1(2).ne.dim2(2)) then
	write(6,*)'Dimensions de les matrius diferents !!!!'
	stop
	endif

	nx=dim1(1)
	ny=dim1(2)
	write(6,*)'Dimensio-matrius... ',nx,' x ',ny	
	write(6,*)' '

	open(unit=1,file=ubicacio1(1:n_ubic1) // char(47) // 
     &arxiu_1(1:n_a1),status='old')
  	open(unit=2,file=ubicacio2(1:n_ubic2) // char(47) // 
     &arxiu_2(1:n_a2),status='old')
   	
	do i=1,nx
	read(1,40)(matriu1(i,j),j=1,ny)
	read(2,40)(matriu2(i,j),j=1,ny)
  40	format(300(f16.6))
	enddo

	close(unit=1)
	close(unit=2)

	resta_min=10.**(10.)
	resta_max=-10.**(10.)
	mig1=0.
	migquad1=0.
	mig2=0.
	migquad2=0.
	mig12=0.
	migquad12=0.
	migresta=0.
	migquadresta=0.

	write(6,*)'Calculant la diferencia...'
	write(6,*)'Matriu 2 - Matriu 1'
	write(6,*)' '
	do i=1,nx
	do j=1,ny
	mig1=mig1+matriu1(i,j)
	migquad1=migquad1+matriu1(i,j)*matriu1(i,j)
	mig2=mig2+matriu2(i,j)
	migquad2=migquad2+matriu2(i,j)*matriu2(i,j)
	mig12=mig12+matriu1(i,j)*matriu2(i,j)
	migquad12=migquad12+matriu1(i,j)*matriu2(i,j)*
     &matriu1(i,j)*matriu2(i,j)
	resta(i,j)=matriu2(i,j)-matriu1(i,j)
	migresta=migresta+resta(i,j)
	migquadresta=migquadresta+resta(i,j)*resta(i,j)
	if (resta(i,j).lt.resta_min) resta_min=resta(i,j)
	if (resta(i,j).gt.resta_max) resta_max=resta(i,j)	
	enddo
	enddo

	migresta=migresta/(nx*ny*1.)
	migquadresta=migquadresta/(nx*ny*1.)
	mig1=mig1/(nx*ny*1.)
	migquad1=migquad1/(nx*ny*1.)
	mig2=mig2/(nx*ny*1.)
	migquad2=migquad2/(nx*ny*1.)
	mig12=mig12/(nx*ny*1.)
	migquad12=migquad12/(nx*ny*1.)
	
	varresta=sqrt(migquadresta-migresta*migresta)
	var1=migquad1-mig1*mig1
	var2=migquad2-mig2*mig2
	var12=mig12-mig1*mig2
	
	rmse=sqrt(migquadresta-migresta*migresta)
	corr=sqrt(var12*var12/(var1*var2))

	mae=0.
	do i=1,nx
	do j=1,ny
	mae=mae+abs(migresta-resta(i,j))
	enddo
	enddo
	mae=mae/(nx*ny*1.)

  	open(unit=3,file=ubicacio_sortida(1:n_ubics) // char(47) // 
     &arxiu_sortida(1:n_as),status='unknown')

	write(6,*)'Escrivint arxiu de sortida.... '
	write(6,*)'Resultats respecte el valor de la matriu 2'
	write(6,*)' '
	write(6,*)ubicacio1(1:n_ubic1) // char(47) // 
     &arxiu_sortida(1:n_as)
	write(6,*)' '
	write(6,*)'******************************'
	write(6,*)'Diferencia minima: ',resta_min
	write(6,*)'Diferencia maxima: ',resta_max
	write(6,*)'Mig-resta: ',migresta
	write(6,*)'Migquad-resta: ',migquadresta
	write(6,*)'Varresta: ',varresta
	write(6,*)'mig1: ',mig1
	write(6,*)'migquad1: ',migquad1
	write(6,*)'var1: ',var1
	write(6,*)'mig2: ',mig2
	write(6,*)'migquad2: ',migquad2
	write(6,*)'var2: ',var2
	write(6,*)'mig12: ',mig12
	write(6,*)'migquad12: ',migquad12
	write(6,*)'var12: ',var12
	write(6,*)'mae: ',mae
	write(6,*)'rmse: ',rmse
	write(6,*)'Correlacio: ',corr
		
	do i=1,nx
	write(3,50)(resta(i,j),j=1,ny)
   50	format(300(f20.10,1x))
	enddo

	close(unit=3)
	
	end
	
ccccccc	SUBROUTINA per a determinar la llargada d'un String
	subroutine llargada_String(String,num_caracters,llargada)
	integer num_caracters,llargada
	character*1000 String
	
	llargada=0
	do i=1,num_caracters
 	if ((ichar(String(i:i))-48).ne.(-16)) llargada=llargada+1
	enddo

c	write(6,*)'*****',String(1:llargada),'****'
	end


ccccccc	SUBROUTINA per a llegir les dades i buscar-ne la seva dimensio
	subroutine dimensions(ar,ubic,dx,dy)
	integer i,j
	integer dx,dy,n_ar,n_u,ll,charlin
	character*50 lin
	character*80 ar
	character*200 ubic
	character*5000 linia 
	dimension charlin(50)
	
	n_ar=0
	do i=1,80
	if ((ichar(ar(i:i))-48).ne.(-16)) n_ar=n_ar+1
	enddo
	
	n_u=0
	do i=1,200
	if ((ichar(ubic(i:i))-48).ne.(-16)) n_u=n_u+1
	enddo

	dx=0
	dy=0
	open(unit=3,file=ubic(1:n_u) // char(47) // ar(1:n_ar)
     &,status='old')
	do i=1,1000
	read(3,5,end=100)lin
	dx=dx+1
	enddo
  100	continue
  	close(unit=3)
	
	open(unit=3,file=ubic(1:n_u) // char(47) // ar(1:n_ar)
     &,status='old')
	
	read(3,70)linia
  70 	format(5000a)

	close(unit=3)

	do i=1,50
	charlin(i)=ichar(linia(i:i))-48
	enddo

	ll=0
	do i=2,49
	if ((charlin(i).eq.-16).and.(charlin(i-1).ne.(-16))) then
	ll=i-1
	goto 150
	endif
	enddo
  150	continue

ccc
c	Cas amb espai entre valors
ccc
c	ll=ll+1

	dy=0
	do i=1,1000
	ilin=ll-1+(i-1)*ll
	if ((ichar(linia(ilin:ilin))-48).ne.(-16)) dy=dy+1
	if ((ichar(linia(ilin:ilin))-48).eq.(-16)) goto 200
	enddo
  200	continue

c 	write(6,*)'La dimensio de la matriu ',dx,' x ',dy

   5	format(a50)
	end

ccccccc Programa per a llegir dades d'un arxiu de text
	subroutine llegir(ubic1,ubic2,arx_1,arx_2,ubic_s,arx_s)
	character*10 linia
	character*2 2lletres
	character*80 arx_1,arx_2,arx_s
	character*200 ubic1,ubic2,ubic_s
	
	open(unit=8,file='inf_dif_matriu.txt',status='old')
	
	read(8,10)linia
	read(8,10)linia
	read(8,20)ubic1
	read(8,10)linia
	read(8,10)linia
	read(8,20)ubic2
	read(8,10)linia
	read(8,10)linia
	read(8,80)arx_1
	read(8,10)linia
	read(8,10)linia
	read(8,80)arx_2
	read(8,10)linia
	read(8,10)linia
	read(8,20)ubic_s
	read(8,10)linia
	read(8,10)linia
	read(8,80)arx_s


 80	format(a80)
 10	format(a10)
 20	format(a200)

	write(6,*)'Tenim:'
	write(6,*)'Dimensions matriu de dades: Dim x ',nx,' Dim y ',ny
	write(6,*)'Ubicacio dels arxius: ',ubic1
	write(6,*)'Ubicacio dels arxius: ',ubic2	
	write(6,*)'Arxiu de la 1a matriu: ',arx_1
	write(6,*)'Arxiu de la 2a matriu: ',arx_2
	write(6,*)'Ubicacio arxius de sortida: ',ubic_s
	write(6,*)'Arxiu de sortida de la matriu: ',arx_s

	close(unit=8)
	end
	
	
	
