diff --git a/src/MNH/diag.f90 b/src/MNH/diag.f90
index 2752a91a5809a7d00e57d7c5f270320b1894a2df..634b5d351f9613ed49c336fb76d3365346853bd8 100644
--- a/src/MNH/diag.f90
+++ b/src/MNH/diag.f90
@@ -2,7 +2,7 @@
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
 !MNH_LIC for details. version 1.
-! $Source$ $Revision$ $Date$
+! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/diag.f90,v $ $Revision: 1.3.2.4.2.4.2.5.2.6.2.3.2.6 $ $Date: 2015/11/26 14:55:02 $
 !-----------------------------------------------------------------
 !     ############
       PROGRAM DIAG
@@ -76,6 +76,7 @@
 !!  04/2016     (G.Delautier) replace print by write in OUTPUT LISTING
 !!  06/2016     (G.Delautier) phasage surfex 8
 !!  09/2016      (JP Pinty) Add LIMA
+!!  10/2016      (C.LAC) add LVISI
 !-------------------------------------------------------------------------------
 !
 !*       0.     DECLARATIONS
@@ -222,7 +223,7 @@ NAMELIST/NAM_DIAG/ CISO, LVAR_RS, LVAR_LS,   &
                    XGRID,NBELEV,XELEV,NBRAD,LQUAD,LFALL,LWBSCS,LWREFL,&
                    XREFLMIN,XREFLVDOPMIN,LSNRT,XSNRMIN,&
                    LLIDAR,CVIEW_LIDAR,XALT_LIDAR,XWVL_LIDAR,&
-                   LISOPR,XISOPR,LISOTH,XISOTH, LHU_FLX, LLIMA_DIAG
+                   LISOPR,XISOPR,LISOTH,XISOTH, LHU_FLX,LVISI, LLIMA_DIAG
 !
 NAMELIST/NAM_DIAG_FILE/ YINIFILE,YINIFILEPGD, YSUFFIX
 NAMELIST/NAM_STO_FILE/ CFILES, NSTART_SUPP
@@ -285,6 +286,7 @@ NCAPE=-1
 LBV_FR=.FALSE.
 LRADAR=.FALSE.
 LBLTOP=.FALSE.
+LVISI=.FALSE.
 LVAR_FRC=.FALSE.
 LCHEMDIAG=.FALSE.
 CAERDIAG='CLIM'
diff --git a/src/MNH/ini_posprofilern.f90 b/src/MNH/ini_posprofilern.f90
index fa869c0763ea0a1b1c085e35f53e8c6209f293d4..0608f5764240d2178c0a67387e59ef07934690c4 100644
--- a/src/MNH/ini_posprofilern.f90
+++ b/src/MNH/ini_posprofilern.f90
@@ -72,7 +72,7 @@ END MODULE MODI_INI_POSPROFILER_n
 !!    MODIFICATIONS
 !!    -------------
 !!     P. Tulet 15/01/2002  
-!!
+!!     C.Lac 10/2016  Add visibility diagnostic
 !! --------------------------------------------------------------------------
 !       
 !*      0. DECLARATIONS
@@ -182,6 +182,8 @@ ALLOCATE(TPROFILER%ZZ    (ISTORE,IKU,NUMBPROFILER))
 ALLOCATE(TPROFILER%TH    (ISTORE,IKU,NUMBPROFILER))
 ALLOCATE(TPROFILER%THV   (ISTORE,IKU,NUMBPROFILER))
 ALLOCATE(TPROFILER%RHOD  (ISTORE,IKU,NUMBPROFILER))
+ALLOCATE(TPROFILER%VISI  (ISTORE,IKU,NUMBPROFILER))
+ALLOCATE(TPROFILER%VISIKUN(ISTORE,IKU,NUMBPROFILER))
 ALLOCATE(TPROFILER%RARE  (ISTORE,IKU,NUMBPROFILER))
 ALLOCATE(TPROFILER%R     (ISTORE,IKU,NUMBPROFILER,KRR))
 ALLOCATE(TPROFILER%SV    (ISTORE,IKU,NUMBPROFILER,KSV))
@@ -223,6 +225,8 @@ TPROFILER%ZZ   = XUNDEF
 TPROFILER%TH   = XUNDEF
 TPROFILER%THV  = XUNDEF
 TPROFILER%RHOD = XUNDEF
+TPROFILER%VISI = XUNDEF
+TPROFILER%VISIKUN = XUNDEF
 TPROFILER%RARE = XUNDEF
 TPROFILER%IWV  = XUNDEF
 TPROFILER%ZTD  = XUNDEF
diff --git a/src/MNH/modd_diag_flag.f90 b/src/MNH/modd_diag_flag.f90
index 66f3c6f2be66567753a171d859332fb5ee633c30..852e3741336e88261602e88bdc0c0e6fc3a8aa91 100644
--- a/src/MNH/modd_diag_flag.f90
+++ b/src/MNH/modd_diag_flag.f90
@@ -5,7 +5,7 @@
 !-----------------------------------------------------------------
 !--------------- special set of characters for RCS information
 !-----------------------------------------------------------------
-! $Source$ $Revision$
+! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/modd_diag_flag.f90,v $ $Revision: 1.2.4.1.2.2.10.2.2.1 $
 ! masdev4_8 modd 2008/06/30 15:13:13
 !-----------------------------------------------------------------
 !     #####################
@@ -37,6 +37,7 @@
 !!       J.-P. Chaboureau 15/04/03  add LRAD_SUBG_COND
 !!       L. Leriche 21/04/06 add aqueous phase chemistry LCHAQDIAG
 !!       D.Ricard 2015 : add LMOIST_ES
+!!     C.Lac 10/2016  Add visibility diagnostic
 !!
 !-------------------------------------------------------------------------------
 !
@@ -91,6 +92,7 @@ INTEGER     :: NCAPE       ! CAPE, DCAPE, CIN, CAPEMAX, CINMAX
 LOGICAL     :: LBV_FR
 LOGICAL     :: LRADAR
 LOGICAL     :: LBLTOP
+LOGICAL     :: LVISI 
 LOGICAL     :: LTRAJ       ! to compute trajectories
 LOGICAL     :: LCHEMDIAG = .FALSE.  ! flag for chemistry
 CHARACTER (LEN=4) :: CAERDIAG  ! aerosols optical thickness type
diff --git a/src/MNH/modd_type_profiler.f90 b/src/MNH/modd_type_profiler.f90
index 4bd4f703bd2ea5e8b04147a674a801dec5cb4b96..9c14fdda428a3484fe34389787de9aef2f143bfe 100644
--- a/src/MNH/modd_type_profiler.f90
+++ b/src/MNH/modd_type_profiler.f90
@@ -5,7 +5,7 @@
 !-----------------------------------------------------------------
 !--------------- special set of characters for RCS information
 !-----------------------------------------------------------------
-! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/modd_type_profiler.f90,v $ $Revision: 1.2.4.1.2.1.10.2.2.2 $
+! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/modd_type_profiler.f90,v $ $Revision: 1.2.4.1.2.1.10.2.2.1 $
 ! MASDEV4_7 modd 2006/06/27 12:27:06
 !-----------------------------------------------------------------
 !     ############################
@@ -36,6 +36,7 @@ IMPLICIT NONE
 !!    MODIFICATIONS
 !!    -------------
 !!      Original    15/01/02
+!!     C.Lac 10/2016  Add visibility diagnostic
 !-------------------------------------------------------------------------------
 !
 !*       0.   DECLARATIONS
@@ -77,6 +78,8 @@ REAL, DIMENSION(:,:,:),   POINTER :: ZZ=>NULL()       ! altitude(n)
 REAL, DIMENSION(:,:,:),   POINTER :: TKE=>NULL()      ! tke(n)
 REAL, DIMENSION(:,:,:),   POINTER :: TH=>NULL()       ! th(n)
 REAL, DIMENSION(:,:,:),   POINTER :: THV=>NULL()      ! thv(n)
+REAL, DIMENSION(:,:,:),   POINTER :: VISI=>NULL()     ! VISI(n)
+REAL, DIMENSION(:,:,:),   POINTER :: VISIKUN=>NULL()  ! VISI KUNKEL(n)
 REAL, DIMENSION(:,:,:),   POINTER :: RARE=>NULL()     ! radar reflectivity (n)
 REAL, DIMENSION(:,:,:),   POINTER :: RHOD=>NULL()     ! density of dry air/moist air
 REAL, DIMENSION(:,:,:,:), POINTER :: R=>NULL()        ! r*(n)
@@ -84,7 +87,7 @@ REAL, DIMENSION(:,:,:,:), POINTER :: SV=>NULL()       ! Sv*(n)
 REAL, DIMENSION(:,:,:,:), POINTER :: AER=>NULL()      ! AER*(n) aerosol extinction
 REAL, DIMENSION(:,:),     POINTER :: DATIME=>NULL()   ! record for diachro
 !
-REAL, DIMENSION(:,:),     POINTER :: T2M=>NULL()      ! 2 m air temperature (°C)
+REAL, DIMENSION(:,:),     POINTER :: T2M=>NULL()      ! 2 m air temperature (°C)
 REAL, DIMENSION(:,:),     POINTER :: Q2M=>NULL()      ! 2 m humidity (kg/kg)
 REAL, DIMENSION(:,:),     POINTER :: HU2M=>NULL()     ! 2 m relative humidity (%)
 REAL, DIMENSION(:,:),     POINTER :: ZON10M=>NULL()   ! 10 m zonal wind (m/s)  
diff --git a/src/MNH/profilern.f90 b/src/MNH/profilern.f90
index 269caeb1cdce0381ef65dae7f32cef4c7f88750b..64ed6ad1b478cb83c2d5643dccaf71e392798736 100644
--- a/src/MNH/profilern.f90
+++ b/src/MNH/profilern.f90
@@ -5,7 +5,7 @@
 !-----------------------------------------------------------------
 !--------------- special set of characters for RCS information
 !-----------------------------------------------------------------
-! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/profilern.f90,v $ $Revision: 1.3.4.1.2.1.10.2.2.3 $
+! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/profilern.f90,v $ $Revision: 1.3.4.1.2.1.10.2.2.4 $
 ! MASDEV4_7 profiler 2006/06/01 09:51:49
 !-----------------------------------------------------------------
 !      ##########################
@@ -18,7 +18,7 @@ INTERFACE
                             TPDTEXP, TPDTMOD, TPDTSEG, TPDTCUR,   &
                             PXHAT, PYHAT, PZ,PRHODREF,            &
                             PU, PV, PW, PTH, PR, PSV, PTKE,       &
-                         PTS,PP, PAER, PCLDFR, PCIT)
+                            PTS,PP, PAER, PCLDFR, PCIT)
 !
 USE MODD_TYPE_DATE
 !
@@ -58,7 +58,7 @@ END MODULE MODI_PROFILER_n
                             TPDTEXP, TPDTMOD, TPDTSEG, TPDTCUR,   &
                             PXHAT, PYHAT, PZ,PRHODREF,            &
                             PU, PV, PW, PTH, PR, PSV, PTKE,       &
-                         PTS,PP, PAER, PCLDFR, PCIT)
+                            PTS, PP, PAER, PCLDFR, PCIT)
 !     ########################################################
 !
 !
@@ -93,6 +93,7 @@ END MODULE MODI_PROFILER_n
 !!     Original 15/02/2002
 !!     March 2013 : C.Lac : Corrections for 1D + new fields (RARE,THV,DD,FF)
 !!     April 2014 : C.Lac : Call RADAR only if ICE3   
+!!     C.Lac 10/2016  Add visibility diagnostic
 !!
 !!
 !! --------------------------------------------------------------------------
@@ -109,6 +110,7 @@ USE MODD_CST
 USE MODD_GRID
 USE MODD_DIAG_IN_RUN
 USE MODD_CONF
+USE MODD_NSV
 !
 USE MODE_ll
 !
@@ -206,9 +208,6 @@ REAL                       :: ZZHD_PROFILER ! ZHD at station location
 REAL                       :: ZZWD_PROFILER ! ZWD at station location
 REAL                       :: ZZHDR         ! ZHD correction at station location
 REAL                       :: ZZWDR         ! ZWD correction at station location
-!-------- Calculation parameters --------
-REAL ::  ZK1,ZK2,ZK3            ! k1, k2 and K3 atmospheric refractivity constants
-REAL  :: ZRDSRV                 ! XRD/XRV
 !
 INTEGER                    :: IINFO_ll    ! return code
 INTEGER                    :: ILUOUT      ! logical unit
@@ -218,6 +217,9 @@ INTEGER                    :: I           ! loop for stations
 REAL,DIMENSION(SIZE(PTH,1),SIZE(PTH,2))              :: ZZTD,ZZHD,ZZWD
 REAL,DIMENSION(SIZE(PTH,1),SIZE(PTH,2),SIZE(PTH,3))  :: ZTEMP,ZRARE,ZTHV,ZTEMPV
 REAL,DIMENSION(SIZE(PTH,1),SIZE(PTH,2),SIZE(PTH,3))  :: ZWORK32,ZWORK33,ZWORK34,ZCIT
+REAL,DIMENSION(SIZE(PTH,1),SIZE(PTH,2),SIZE(PTH,3))  :: ZVISI,ZVISIKUN
+REAL ::  ZK1,ZK2,ZK3            ! k1, k2 and K3 atmospheric refractivity constants
+REAL  :: ZRDSRV                 ! XRD/XRV
 !----------------------------------------------------------------------------
 !
 !*      2.   PRELIMINARIES
@@ -404,6 +406,19 @@ ZTHV(:,:,:) = PTH(:,:,:) / (1.+WATER_SUM(PR(:,:,:,:)))*(1.+PR(:,:,:,1)/ZRDSRV)
 ! virtual temperature
 ZTEMPV(:,:,:)=ZTHV(:,:,:)*(PP(:,:,:)/ XP00) **(XRD/XCPD)
 CALL GPS_ZENITH_GRID(PR(:,:,:,1),ZTEMP,PP,ZZTD,ZZHD,ZZWD)
+! Kunkel formulation
+IF (SIZE(PR,2) >= 2) THEN
+  WHERE ( PR(:,:,:,2) /=0 )
+    ZVISIKUN(:,:,:) =0.027/(PR(:,:,:,2)*PRHODREF(:,:,:))**0.88
+  END WHERE
+END IF
+! Gultepe formulation
+IF ((SIZE(PR,2) >= 2) .AND. NSV_C2R2END /= 0 ) THEN 
+  WHERE ( (PR(:,:,:,2) /=0. ) .AND. (PSV(:,:,:,NSV_C2R2BEG+1) /=0. ) )
+    ZVISI(:,:,:) =1.002/(PR(:,:,:,2)*PRHODREF(:,:,:)*PSV(:,:,:,NSV_C2R2BEG+1))**0.6473
+  END WHERE
+END IF
+!
 IF (GSTORE) THEN
  IF (TPROFILER%TIME(IN) /= XUNDEF) THEN
   DO I=1,NUMBPROFILER
@@ -493,6 +508,8 @@ IF (GSTORE) THEN
       TPROFILER%W   (IN,:,I) = PROFILER_INTERP(PW)
       TPROFILER%TH  (IN,:,I) = PROFILER_INTERP(PTH)
       TPROFILER%THV (IN,:,I) = PROFILER_INTERP(ZTHV)
+      TPROFILER%VISI(IN,:,I) = PROFILER_INTERP(ZVISI)
+      TPROFILER%VISIKUN(IN,:,I) = PROFILER_INTERP(ZVISIKUN)
       TPROFILER%ZZ  (IN,:,I) = ZZ(:)
       TPROFILER%RHOD(IN,:,I) = ZRHOD(:)
       IF (SIZE(PR,4) == 6) TPROFILER%RARE(IN,:,I) = PROFILER_INTERP(ZRARE)
@@ -574,6 +591,8 @@ IF (GSTORE) THEN
   CALL DISTRIBUTE_PROFILER(TPROFILER%ZZ  (IN,JK,I))
   CALL DISTRIBUTE_PROFILER(TPROFILER%TH  (IN,JK,I))
   CALL DISTRIBUTE_PROFILER(TPROFILER%THV (IN,JK,I))
+  CALL DISTRIBUTE_PROFILER(TPROFILER%VISI(IN,JK,I))
+  CALL DISTRIBUTE_PROFILER(TPROFILER%VISIKUN(IN,JK,I))
   CALL DISTRIBUTE_PROFILER(TPROFILER%RHOD(IN,JK,I))
   CALL DISTRIBUTE_PROFILER(TPROFILER%RARE(IN,JK,I))
   CALL DISTRIBUTE_PROFILER(TPROFILER%IWV(IN,I))
diff --git a/src/MNH/write_profilern.f90 b/src/MNH/write_profilern.f90
index 8f48f7b68924b69bedc215a32005c326195b52ff..b2b24c93f071aeac2b082e0d9610e2fb0add8ba3 100644
--- a/src/MNH/write_profilern.f90
+++ b/src/MNH/write_profilern.f90
@@ -60,6 +60,7 @@ END MODULE MODI_WRITE_PROFILER_n
 !!    -------------
 !!     Original 15/02/2002
 !!     2016 : G.DELAUTIER : LIMA
+!!              Oct, 2016 (C.Lac) Add visibility diagnostics for fog
 !!
 !! --------------------------------------------------------------------------
 !       
@@ -153,7 +154,7 @@ IF (TPROFILER%X(II)==XUNDEF) RETURN
 IF (TPROFILER%Y(II)==XUNDEF) RETURN
 IKU = SIZE(TPROFILER%W,2)    !nbre de niveaux sur la verticale SIZE(TPROFILER%W,2)
 !
-IPROC = 22 + SIZE(TPROFILER%R,4) + SIZE(TPROFILER%SV,4)
+IPROC = 24 + SIZE(TPROFILER%R,4) + SIZE(TPROFILER%SV,4)
 IF (LDIAG_IN_RUN) IPROC = IPROC + 13
 IF (LORILAM) IPROC = IPROC + JPMODE*3
 IF (LDUST) IPROC = IPROC + NMODE_DST*3
@@ -190,6 +191,18 @@ YCOMMENT (JPROC) = 'Virtual Potential temperature'
 ZWORK6 (1,1,IK,:,1,JPROC) = TPROFILER%THV(:,IK,II)
 !
 JPROC = JPROC + 1
+YTITLE   (JPROC) = 'VISI'
+YUNIT    (JPROC) = 'km'
+YCOMMENT (JPROC) = 'Visibility'
+ZWORK6 (1,1,IK,:,1,JPROC) = TPROFILER%VISI(:,IK,II)
+!
+JPROC = JPROC + 1
+YTITLE   (JPROC) = 'VISIKUN'
+YUNIT    (JPROC) = 'km'
+YCOMMENT (JPROC) = 'Visibility Kunkel'
+ZWORK6 (1,1,IK,:,1,JPROC) = TPROFILER%VISIKUN(:,IK,II)
+!
+JPROC = JPROC + 1
 YTITLE   (JPROC) = 'RARE'
 YUNIT    (JPROC) = 'dBz'
 YCOMMENT (JPROC) = 'Radar reflectivity'