c**********************************************************************

	SUBROUTINE lum_teff_VIF
	
c routine PUBLIC du mod_photo

c determination des luminosites, temperatures effectives et ecarts
c a partir de pi, V, I, [Fe/H]

c Auteur: P. Morel, Laboratoire Lagrange, O.C.A., CESAM2k

c------------------------------------------------------------------

	USE mod_kind
		
	IMPLICIT NONE
	
	REAL (kind=dp), DIMENSION(2,4) :: jac
	REAL (kind=dp), DIMENSION(4) ::	x, dx	    	
	REAL (kind=dp), DIMENSION(2) :: f, f0, sig
	REAL (kind=dp),	PARAMETER :: h=1.d-5, unph=1.d0+h
	REAL (kind=dp) :: v, sv, i, si, para, spara, fesh, sfesh, loge,
	1 va, ia, stor, stor0, dstor
	
	INTEGER	:: k,j
	
c----------------------------------------------------------------------
	
2000	FORMAT(8es10.3)

 	loge=log10(EXP(1.d0))
 	
	PRINT*
	PRINT*,'calcul de L a partir de para, V, I, [Fe/H]'
	PRINT*
	
	B1: DO	
	 PRINT*,'entrer pi, sigma pi (mas)Ex: 752.d0 2.7d0, sortie -1 0'
	 READ*,para,spara
	 if(para <= 0.d0)EXIT B1
	 PRINT*
	 WRITE(*,"('para=',es10.3,' +/-',es10.3)")para,spara

	 PRINT*
	 PRINT*,'entrer magnitude V, sigma V, Exemple: -0.01d0 0.005d0'
	 READ*,v,sv
	 PRINT*
	 WRITE(*,"('V=',es10.3,' +/-',es10.3)")v,sv

	 PRINT*
	 PRINT*,'entrer magnitude I, sigma I, Exemple: 0.69d0 0.01d0'
	 READ*,i,si 
	 PRINT*	
	 WRITE(*,"('I=',es10.3,' +/-',es10.3)")i,si
	
	 PRINT*	
	 PRINT*,'entrer [Fe/H], sigma [Fe/H], Exemple: 0.2d0 0.02d0'
	 READ*,fesh,sfesh 		
	 PRINT*	
	 WRITE(*,"('[Fe/H=]',es10.3,' +/-',es10.3)")fesh,sfesh

c variables muettes pour le calcul du jacobien
	 x(1)=para   ; x(2)=v   ; x(3)=i  ; x(4)=fesh
	 dx(1)=spara ; dx(2)=sv ;dx(3)=si ; dx(4)=sfesh

c magnitudes absolues
	 va=m_abs(x(2),x(1)) ; ia=m_abs(x(3),x(1)) 

	 PRINT*
	 WRITE(*,"('Magnitudes absolues')")	
	 WRITE(*,"('I=',es10.3,', V=',es10.3)")ia,va
	 PRINT*
			
	 CALL tl_vif(va,ia,x(4),f0(1),f0(2))
					
	 PRINT*
c	 WRITE(*,2000)f0(1),f0(2)
	 WRITE(*,"('Teff=',es10.3,', L/Lsol=',es10.3)")(f0(k),k=1,2)
	 PRINT*
c	 pause
		
	 DO k=1,4
	  stor0=x(k)
	  stor=stor0*unph
	  if(stor == 0.d0)stor=h
	  dstor=stor-stor0
	  x(k)=stor	
	  va=m_abs(x(2),x(1)) ; ia=m_abs(x(3),x(1))
	  CALL tl_vif(va,ia,x(4),f(1),f(2))
	  x(k)=stor0
	  jac(:,k)=(f-f0)/dstor
	 ENDDO	 
	 PRINT* ; PRINT*,'jacobien' ; WRITE(*,"(4es10.3)")(jac(:,j),j=1,4)
	 	
	 sig=MATMUL(ABS(jac),dx)
	 
c en Log	
	 PRINT*
	 PRINT*,'en Log'
	 WRITE(*,1)LOG10(f0(1)),loge*sig(1)/f0(1)
1	 FORMAT(1x,'Log Teff =',es10.3,' +/-',es10.3)
	 WRITE(*,3)LOG10(f0(2)),loge*sig(2)/f0(2)
3	 FORMAT('Log L/Lsol =',es10.3,' +/-',es10.3)

c en valeurs naturelles
	 PRINT*
	 PRINT*,'en valeurs naturelles'
	 WRITE(*,5)f0(1),sig(1)
5	 FORMAT('Teff =',es10.3,' +/-',es10.3)
	 WRITE(*,7)f0(2),sig(2)
7	 FORMAT('L/Lsol =',es10.3,' +/-',es10.3,/)
	ENDDO b1

	STOP
	
	END SUBROUTINE lum_teff_VIF