zamaniM
Programmer
- Jul 10, 2011
- 12
Hi. i have a problem,i've wrote a cod for calculation K3_interpolated based on tetaS_radian_interpolated .in fact at first code calculates tetaS_radian for 24 hours of day (for 14 days) and then interpolates tetaS_radian for each 1/60 min ( 86400=24*60*1/(1/60))).
the problem is about tetaS_radian_interpolated, the results of tetaS_radian_interpolated are zero and i don't know what's the reason.for check the problem, i've replaced equations of tetaS_radian_interpolated calculation by a constant. Even with this situation results are zero.
do you have any idea?
PROGRAM CalculationOfk_3
IMPLICIT NONE
!************ Variabels ********************
INTEGER ::j,i,D_l,ii
INTEGER,DIMENSION(14) :_J
REAL ::fi_radian
REAL,PARAMETER :i=3.141592654,year=2010
REAL,DIMENSION(25,14)::tetaS,tetaS_radian
REAL,DIMENSION(86400,2:15)::K3_interpolated
REAL,DIMENSION(25) ::time_hour,t_s
REAL,DIMENSION(25) ::H_a,Ha_radian
REAL,DIMENSION(14) ::N_J_D,EpsilonOB_radian,Delta_radian,landaEC
REAL,DIMENSION(14) ::landaEC_radian,g_M_radian,L_M,L_M_radian
REAL,DIMENSION(86400,15)::tetaS_radian_interpolated
REAL ,DIMENSION (24,14)::slop_tetaS_radian
REAL::time,dt,q,qqq
!***************variabels*******************************
dt=1/60
!******************************* OPEN ******
OPEN(1,file='D_J.txt',status='old')
OPEN(2,file='tetaS_radian_interpolated .txt')
OPEN(9,file='k3_interpolated.txt')
OPEN(11,file='tetaS.txt')
OPEN(12,file='tetaS_radian.txt')
!**************************************
!******************************** READ ****************
DO j=1,14
READ(1,*)(D_J(j))
END DO
!*************************************************
!********** Calulation of Photolysis Rate ************
fi_radian=(35.0+(48.0/60.0)+(0.214/3600.0))*(pi/180)
DO j=1,14
DO i=1,25
time_hour(i)=i-1
t_s(i)=(time_hour(i)-12)*3600
H_a(i)=((2*pi*t_s(i))/86400)*(180/pi)
IF(year<2001)THEN
D_L=INT((year-2001)/4)-1
ELSE
D_L=INT((year-2001)/4)
END IF
N_J_D(j)=364.5+((year-2001)*365)+D_L+D_J(j)
EpsilonOB_radian(j)=(23.49-(0.0000004*N_J_D(j)))*(pi/180)
L_M(j)=280.460+(0.9856474*N_J_D(j))
L_M_radian(j)=(L_M(j))*(pi/180)
g_M_radian(j)=(357.528+(0.9856003*N_J_D(j)))*(pi/180)
q=(0.020*SIN(2*g_M_radian(j)))
landaEC(j)=L_M(j)+(1.915*SIN(g_M_radian(j)))+q
landaEC_radian(j)=( landaEC(j))*(pi/180)
Ha_radian(i)=H_a(i)*(pi/180)
q=SIN(landaEC_radian(j))
Delta_radian(j)=(ASIN(SIN(EpsilonOB_radian(j))*q))
q=(COS(fi_radian)*COS(Delta_radian(j))*COS(Ha_radian(i)))
tetaS_radian(i,j)=ACOS((SIN(fi_radian)*SIN(Delta_radian(j)))+q)
tetaS(i,j)=tetaS_radian(i,j)*(180/pi)
END DO
END DO
DO j=1,14
DO i=1,24
q=(tetaS_radian(i+1,j)-tetaS_radian(i,j))
slop_tetaS_radian(i,j)=q/60
END DO
END DO
DO j=1,14
DO i=1,24
DO ii=1,60/dt
time=((ii-1)*dt)+((i-1)*60)
q=slop_tetaS_radian(i,j)*(time-((i-1)*60))
qqq=tetaS_radian(i,j)
tetaS_radian_interpolated(((i-1)*60/dt)+ii,j+1)=q+qqq
tetaS_radian_interpolated(((i-1)*60/dt)+ii,1)=1
END DO
END DO
END DO
DO j=2,15
DO i=1,86400
k3_interpolated(i,j)=EXP(-0.575/SIN(tetaS_radian_interpolated (i,j)))
END DO
END DO
DO i=1,86400
WRITE(2,"(15(1x,1pe15.7))") (tetaS_radian_interpolated (i,j),j=1,15)
WRITE(9,"(14(1x,1pe14.7))")(k3_interpolated(i,j),j=2,15)
END DO
!********************************
!******** WRITE **********************
DO i=1,25
WRITE(11,"(14(1x,1pe14.7))")(tetaS(i,j),j=1,14)
WRITE(12,"(14(1x,1pe14.7))")(tetaS_radian(i,j),j=1,14)
END DO
!**************************************
END PROGRAM CalculationOfk_3
the problem is about tetaS_radian_interpolated, the results of tetaS_radian_interpolated are zero and i don't know what's the reason.for check the problem, i've replaced equations of tetaS_radian_interpolated calculation by a constant. Even with this situation results are zero.
do you have any idea?
PROGRAM CalculationOfk_3
IMPLICIT NONE
!************ Variabels ********************
INTEGER ::j,i,D_l,ii
INTEGER,DIMENSION(14) :_J
REAL ::fi_radian
REAL,PARAMETER :i=3.141592654,year=2010
REAL,DIMENSION(25,14)::tetaS,tetaS_radian
REAL,DIMENSION(86400,2:15)::K3_interpolated
REAL,DIMENSION(25) ::time_hour,t_s
REAL,DIMENSION(25) ::H_a,Ha_radian
REAL,DIMENSION(14) ::N_J_D,EpsilonOB_radian,Delta_radian,landaEC
REAL,DIMENSION(14) ::landaEC_radian,g_M_radian,L_M,L_M_radian
REAL,DIMENSION(86400,15)::tetaS_radian_interpolated
REAL ,DIMENSION (24,14)::slop_tetaS_radian
REAL::time,dt,q,qqq
!***************variabels*******************************
dt=1/60
!******************************* OPEN ******
OPEN(1,file='D_J.txt',status='old')
OPEN(2,file='tetaS_radian_interpolated .txt')
OPEN(9,file='k3_interpolated.txt')
OPEN(11,file='tetaS.txt')
OPEN(12,file='tetaS_radian.txt')
!**************************************
!******************************** READ ****************
DO j=1,14
READ(1,*)(D_J(j))
END DO
!*************************************************
!********** Calulation of Photolysis Rate ************
fi_radian=(35.0+(48.0/60.0)+(0.214/3600.0))*(pi/180)
DO j=1,14
DO i=1,25
time_hour(i)=i-1
t_s(i)=(time_hour(i)-12)*3600
H_a(i)=((2*pi*t_s(i))/86400)*(180/pi)
IF(year<2001)THEN
D_L=INT((year-2001)/4)-1
ELSE
D_L=INT((year-2001)/4)
END IF
N_J_D(j)=364.5+((year-2001)*365)+D_L+D_J(j)
EpsilonOB_radian(j)=(23.49-(0.0000004*N_J_D(j)))*(pi/180)
L_M(j)=280.460+(0.9856474*N_J_D(j))
L_M_radian(j)=(L_M(j))*(pi/180)
g_M_radian(j)=(357.528+(0.9856003*N_J_D(j)))*(pi/180)
q=(0.020*SIN(2*g_M_radian(j)))
landaEC(j)=L_M(j)+(1.915*SIN(g_M_radian(j)))+q
landaEC_radian(j)=( landaEC(j))*(pi/180)
Ha_radian(i)=H_a(i)*(pi/180)
q=SIN(landaEC_radian(j))
Delta_radian(j)=(ASIN(SIN(EpsilonOB_radian(j))*q))
q=(COS(fi_radian)*COS(Delta_radian(j))*COS(Ha_radian(i)))
tetaS_radian(i,j)=ACOS((SIN(fi_radian)*SIN(Delta_radian(j)))+q)
tetaS(i,j)=tetaS_radian(i,j)*(180/pi)
END DO
END DO
DO j=1,14
DO i=1,24
q=(tetaS_radian(i+1,j)-tetaS_radian(i,j))
slop_tetaS_radian(i,j)=q/60
END DO
END DO
DO j=1,14
DO i=1,24
DO ii=1,60/dt
time=((ii-1)*dt)+((i-1)*60)
q=slop_tetaS_radian(i,j)*(time-((i-1)*60))
qqq=tetaS_radian(i,j)
tetaS_radian_interpolated(((i-1)*60/dt)+ii,j+1)=q+qqq
tetaS_radian_interpolated(((i-1)*60/dt)+ii,1)=1
END DO
END DO
END DO
DO j=2,15
DO i=1,86400
k3_interpolated(i,j)=EXP(-0.575/SIN(tetaS_radian_interpolated (i,j)))
END DO
END DO
DO i=1,86400
WRITE(2,"(15(1x,1pe15.7))") (tetaS_radian_interpolated (i,j),j=1,15)
WRITE(9,"(14(1x,1pe14.7))")(k3_interpolated(i,j),j=2,15)
END DO
!********************************
!******** WRITE **********************
DO i=1,25
WRITE(11,"(14(1x,1pe14.7))")(tetaS(i,j),j=1,14)
WRITE(12,"(14(1x,1pe14.7))")(tetaS_radian(i,j),j=1,14)
END DO
!**************************************
END PROGRAM CalculationOfk_3