875c875
<       use chapman_module,only: chi     ! was t2 in old sub opflux
---
>       use chapman_module,only: chi,slt     ! was t2 in old sub opflux
885,886c885,886
<      |  phid =  2.0e8,
<      |  phin = -2.0e8,
---
>      |  phid =  2.0e8, !To balance the reduce value in night *2.0e8,
>      |  phin = -1.0e8,
904a905,916
> !          if (chi(i,lat)-0.5*pi >= 0.) then
> !            opflux(i,lat) = fen(i)
> !          else
> !            opflux(i,lat) = fed(i)
> !          endif
> !          if ((chi(i,lat)*rtd-80.)*(chi(i,lat)*rtd-100.) < 0.) 
> !     |      opflux(i,lat) = .5*(fed(i)+fen(i))+.5*(fed(i)-fen(i))*
> !     |        cos(pi*(chi(i,lat)*rtd-80.)/20.)
> !
> !Revised new local time variation function to match the F3/C profiles
> !Added by ITL, 2012.05.19 at HAO
> !Only reduces the flux during early morning!   
906c918,922
<             opflux(i,lat) = fen(i)
---
>             if (slt(i,lat) < 12.) then
>               opflux(i,lat) = 0.25*fen(i)
>             else
>               opflux(i,lat) = fen(i)
>             endif
910,913c926,951
<           if ((chi(i,lat)*rtd-80.)*(chi(i,lat)*rtd-100.) < 0.) 
<      |      opflux(i,lat) = .5*(fed(i)+fen(i))+.5*(fed(i)-fen(i))*
<      |        cos(pi*(chi(i,lat)*rtd-80.)/20.)
< !
---
> 
> !Day-night terminator zone 
>           if ((chi(i,lat)*rtd-80.)*(chi(i,lat)*rtd-100.) < 0.) then
>             if (slt(i,lat) < 12.) then
>               opflux(i,lat) = .5*(fed(i)+0.5*fen(i))+
>      |                        .5*(fed(i)-0.5*fen(i))*
>      |                        cos(pi*(chi(i,lat)*rtd-80.)/20.)
>             else
>               opflux(i,lat) = .5*(fed(i)+fen(i))+.5*(fed(i)-fen(i))*
>      |                        cos(pi*(chi(i,lat)*rtd-80.)/20.)
>             endif
>           endif
> 
> !midnight zone
>           if ((chi(i,lat)*rtd-160.)*(chi(i,lat)*rtd-180.) <= 0.) then
>             if (slt(i,lat) < 12.) then
>               opflux(i,lat) = 0.75*fen(i)+.25*fen(i)*
>      |                        cos(0.5*pi+pi*(chi(i,lat)*rtd-160.)/20.)
>             else
>               opflux(i,lat) = .5*(fen(i)+0.75*fen(i))+
>      |                        .5*(fen(i)-0.75*fen(i))*
>      |                        cos(pi*(chi(i,lat)*rtd-160.)/20.)
>             endif
>           endif 
> 
>