Index: /LMDZ5/trunk/DefLists/context_lmdz.xml
===================================================================
--- /LMDZ5/trunk/DefLists/context_lmdz.xml	(revision 2954)
+++ /LMDZ5/trunk/DefLists/context_lmdz.xml	(revision 2955)
@@ -127,4 +127,10 @@
         <axis id="ReffLiq" />
      </grid>
+     <grid id="grid4Dcol2">
+        <domain id="dom_glo" />
+        <axis id="height" />
+        <axis id="column" />
+     </grid>
+
   </grid_definition>
 
Index: /LMDZ5/trunk/DefLists/cosp_output_nl.txt
===================================================================
--- /LMDZ5/trunk/DefLists/cosp_output_nl.txt	(revision 2954)
+++ /LMDZ5/trunk/DefLists/cosp_output_nl.txt	(revision 2955)
@@ -69,4 +69,15 @@
   Lclmcalipsoun=.true.,
   Lcltcalipsoun=.true.,
+  ! OPAQ variables
+  Lclopaquecalipso=.true.,
+  Lclthincalipso=.true.,
+  Lclzopaquecalipso=.true.,
+  Lclcalipsoopaque=.true.,
+  Lclcalipsothin=.true.,
+  Lclcalipsozopaque=.true.,
+  Lclcalipsoopacity=.true.,
+  ! TIBO variables
+  Lproftemp=.true.,
+  LprofSR=.true.,
   !- ISCCP
   Lalbisccp=.true.,
Index: /LMDZ5/trunk/DefLists/field_def_lmdz.xml
===================================================================
--- /LMDZ5/trunk/DefLists/field_def_lmdz.xml	(revision 2954)
+++ /LMDZ5/trunk/DefLists/field_def_lmdz.xml	(revision 2955)
@@ -796,7 +796,16 @@
       <field id="clcalipsotmpice" long_name="Lidar Ice-Phase Cloud Fraction (532 nm)" unit="1" axis_ref="temp"  detect_missing_value=".true." />
       <field id="clcalipsotmpun"  long_name="Undefined-Phase Total Cloud Fraction"   unit="1" axis_ref="temp"  detect_missing_value=".true." />
+      <field id="clopaquecalipso"  long_name="Lidar Opaque Cloud Fraction"   unit="1"  detect_missing_value=".true." />
+      <field id="clthincalipso"    long_name="Lidar Thin Cloud Fraction"   unit="1"  detect_missing_value=".true." />
+      <field id="clzopaquecalipso" long_name="Lidar mean altitude of opacity"   unit="m"  detect_missing_value=".true." />
+      <field id="clcalipsoopaque"  long_name="Lidar Opaque Cloud profile Fraction (532 nm)" unit="1" axis_ref="height"  detect_missing_value=".true." />
+      <field id="clcalipsothin"    long_name="Lidar Thin Cloud profile Fraction (532 nm)" unit="1" axis_ref="height"  detect_missing_value=".true." />
+      <field id="clcalipsozopaque" long_name="Lidar z_opaque Fraction (532 nm)" unit="1" axis_ref="height"  detect_missing_value=".true." />
+      <field id="clcalipsoopacity" long_name="Lidar opacity Fraction (532 nm)" unit="1" axis_ref="height"  detect_missing_value=".true." />
+      <field id="proftemp"         long_name="Temperature profile (40 lev)" unit="K" axis_ref="height"  detect_missing_value=".true." />
     </field_group>  
       <field id="cfad_lidarsr532" long_name="Lidar Scattering Ratio CFAD (532 nm)" unit="1" freq_op="3h" grid_ref="grid4Dsrbin"  detect_missing_value=".true." />
       <field id="atb532" long_name="Lidar Attenuated Total Backscatter (532 nm)" unit="1" freq_op="3h" grid_ref="grid4Dcol"  detect_missing_value=".true." />
+      <field id="profSR" long_name="Scattering Ratio profiles" unit="1" freq_op="3h" grid_ref="grid4Dcol2"  detect_missing_value=".true." />
 
     <field_group id="fields_COSP_CLOUDSAT" domain_ref="dom_glo" freq_op="3h">
Index: /LMDZ5/trunk/DefLists/file_def_histdayCOSP_lmdz.xml
===================================================================
--- /LMDZ5/trunk/DefLists/file_def_histdayCOSP_lmdz.xml	(revision 2954)
+++ /LMDZ5/trunk/DefLists/file_def_histdayCOSP_lmdz.xml	(revision 2955)
@@ -34,4 +34,13 @@
                 <field field_ref="atb532" level="10" />
                 <field field_ref="beta_mol532" level="10" />
+                <field field_ref="clopaquecalipso"         level="1" />
+                <field field_ref="clthincalipso"           level="1" />
+                <field field_ref="clzopaquecalipso"        level="1" />
+                <field field_ref="clcalipsoopaque"   level="1" axis_ref="height" />
+                <field field_ref="clcalipsothin"     level="1" axis_ref="height" />
+                <field field_ref="clcalipsozopaque"  level="1" axis_ref="height" />
+                <field field_ref="clcalipsoopacity"  level="1" axis_ref="height" />
+                <field field_ref="proftemp"  level="1" axis_ref="height"  />
+                <field field_ref="profSR"  level="1" />
               </field_group>
              <!-- VARS CLOUDSAT et/ou Calipso -->
Index: /LMDZ5/trunk/DefLists/file_def_histhfCOSP_lmdz.xml
===================================================================
--- /LMDZ5/trunk/DefLists/file_def_histhfCOSP_lmdz.xml	(revision 2954)
+++ /LMDZ5/trunk/DefLists/file_def_histhfCOSP_lmdz.xml	(revision 2955)
@@ -34,4 +34,13 @@
                 <field field_ref="atb532" level="10" />
                 <field field_ref="beta_mol532" level="10" />
+                <field field_ref="clopaquecalipso"       level="10" />
+                <field field_ref="clthincalipso"         level="10" />
+                <field field_ref="clzopaquecalipso"      level="10" />
+                <field field_ref="clcalipsoopaque"   level="10" axis_ref="height" />
+                <field field_ref="clcalipsothin"     level="10" axis_ref="height" />
+                <field field_ref="clcalipsozopaque"  level="10" axis_ref="height" />
+                <field field_ref="clcalipsoopacity"  level="10" axis_ref="height" />
+                <field field_ref="proftemp"  level="1" axis_ref="height" />
+                <field field_ref="profSR"  level="1" />
               </field_group>
              <!-- VARS CLOUDSAT et/ou Calipso -->
Index: /LMDZ5/trunk/DefLists/file_def_histmthCOSP_lmdz.xml
===================================================================
--- /LMDZ5/trunk/DefLists/file_def_histmthCOSP_lmdz.xml	(revision 2954)
+++ /LMDZ5/trunk/DefLists/file_def_histmthCOSP_lmdz.xml	(revision 2955)
@@ -35,4 +35,13 @@
                 <field field_ref="atb532" level="10" />
                 <field field_ref="beta_mol532" level="10" />
+                <field field_ref="clopaquecalipso"       level="1" />
+                <field field_ref="clthincalipso"         level="1" />
+                <field field_ref="clzopaquecalipso"      level="1" />
+                <field field_ref="clcalipsoopaque"   level="1" axis_ref="height" />
+                <field field_ref="clcalipsothin"     level="1" axis_ref="height" />
+                <field field_ref="clcalipsozopaque"  level="1" axis_ref="height" />
+                <field field_ref="clcalipsoopacity"  level="1" axis_ref="height" />
+                <field field_ref="proftemp"  level="1" axis_ref="height" />
+                <field field_ref="profSR"  level="1" />
               </field_group>
              <!-- VARS CLOUDSAT et/ou Calipso -->
Index: /LMDZ5/trunk/libf/phylmd/cosp/cosp_constants.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/cosp_constants.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/cosp_constants.F90	(revision 2955)
@@ -54,8 +54,8 @@
 
     ! Number of possible output variables
-    integer,parameter :: N_OUT_LIST = 65
-    integer,parameter :: N3D = 10
-    integer,parameter :: N2D = 14
-    integer,parameter :: N1D = 40
+    integer,parameter :: N_OUT_LIST = 74 !OPAQ 65+7 !TIBO 72+2
+    integer,parameter :: N3D = 11                   !TIBO 10+1
+    integer,parameter :: N2D = 19        !OPAQ 14+4 !TIBO 18+1
+    integer,parameter :: N1D = 43        !OPAQ 40+3
 
     ! Value for forward model result from a level that is under the ground
@@ -96,4 +96,5 @@
     ! Other constants
     integer,parameter :: LIDAR_NCAT    =   4
+    integer,parameter :: LIDAR_NTYPE   =   3 !OPAQ
     integer,parameter :: PARASOL_NREFL =   5 ! parasol
     real,parameter,dimension(PARASOL_NREFL) :: PARASOL_SZA = (/0.0, 20.0, 40.0, 60.0, 80.0/)
Index: /LMDZ5/trunk/libf/phylmd/cosp/cosp_output_mod.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/cosp_output_mod.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/cosp_output_mod.F90	(revision 2955)
@@ -103,4 +103,24 @@
          "clcalipsotmp", "Lidar Cloud Fraction", "%", (/ ('', i=1, 3) /))
 
+  TYPE(ctrl_outcosp), SAVE :: o_clopaquecalipso = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &  !OPAQ
+         "clopaquecalipso", "Lidar Opaque Cloud Fraction", "%", (/ ('', i=1, 3) /))             !OPAQ
+  TYPE(ctrl_outcosp), SAVE :: o_clthincalipso = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &    !OPAQ
+         "clthincalipso", "Lidar Thin Cloud Fraction", "%", (/ ('', i=1, 3) /))                 !OPAQ
+  TYPE(ctrl_outcosp), SAVE :: o_clzopaquecalipso = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), & !OPAQ
+         "clzopaquecalipso", "Lidar mean opacity altitude", "m", (/ ('', i=1, 3) /))            !OPAQ
+  TYPE(ctrl_outcosp), SAVE :: o_clcalipsoopaque = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &  !OPAQ
+         "clcalipsoopaque", "Lidar Opaque profile Cloud Fraction", "%", (/ ('', i=1, 3) /))     !OPAQ
+  TYPE(ctrl_outcosp), SAVE :: o_clcalipsothin = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &    !OPAQ
+         "clcalipsothin", "Lidar Thin profile Cloud Fraction", "%", (/ ('', i=1, 3) /))         !OPAQ
+  TYPE(ctrl_outcosp), SAVE :: o_clcalipsozopaque = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), & !OPAQ
+         "clcalipsozopaque", "Lidar z_opaque Fraction", "%", (/ ('', i=1, 3) /))	        !OPAQ
+  TYPE(ctrl_outcosp), SAVE :: o_clcalipsoopacity = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), & !OPAQ
+         "clcalipsoopacity", "Lidar opacity Fraction", "%", (/ ('', i=1, 3) /))	                !OPAQ
+
+  TYPE(ctrl_outcosp), SAVE :: o_proftemp = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &         !TIBO
+         "proftemp", "Temperature profiles (40 lev)", "K", (/ ('', i=1, 3) /))                  !TIBO
+  TYPE(ctrl_outcosp), SAVE :: o_profSR = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &           !TIBO
+         "profSR", "Lidar Scattering Ratio profiles (532 nm)", "1", (/ ('', i=1, 3) /))         !TIBO
+
 ! Radar Cloudsat
   TYPE(ctrl_outcosp), SAVE :: o_cfadDbze94 = ctrl_outcosp((/ .TRUE., .TRUE., .TRUE. /), &
@@ -376,5 +396,5 @@
       CALL histvert(cosp_nidfiles(iff),"pressure2","pressure","mb",7,ISCCP_PC,nvertisccp(iff),"down")
 
-      CALL histvert(cosp_nidfiles(iff),"column","column","count",Ncolumns,column_ax(1:Ncolumns),nvertcol(iff))
+      CALL histvert(cosp_nidfiles(iff),"column","column","count",Ncolumns,column_ax,nvertcol(iff)) !DBUG
 
       CALL histvert(cosp_nidfiles(iff),"temp","temperature","C",LIDAR_NTEMP,LIDAR_PHASE_TEMP,nverttemp(iff))
Index: /LMDZ5/trunk/libf/phylmd/cosp/cosp_output_write_mod.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/cosp_output_write_mod.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/cosp_output_write_mod.F90	(revision 2955)
@@ -104,4 +104,7 @@
       stlidar%lidarcld(ip,k)=missing_val
      endif
+     if(stlidar%proftemp(ip,k).eq.R_UNDEF)then !TIBO
+      stlidar%proftemp(ip,k)=missing_val       !TIBO
+     endif                                     !TIBO
      enddo
 
@@ -113,5 +116,15 @@
       enddo
      enddo
-   enddo
+
+     do ii= 1,Ncolumns                               !TIBO
+      do ip = 1,Npoints                              !TIBO
+!       if(stlidar%profSR(ip,ii,k).eq.R_UNDEF)then    !TIBO
+!        stlidar%profSR(ip,ii,k)=missing_val          !TIBO
+       if(stlidar%profSR(ip,k,ii).eq.R_UNDEF)then    !TIBO2
+        stlidar%profSR(ip,k,ii)=missing_val          !TIBO2
+       endif                                         !TIBO
+      enddo                                          !TIBO
+     enddo                                           !TIBO
+  enddo
 
   do ip = 1,Npoints
@@ -138,4 +151,12 @@
    enddo
 
+   do k = 1,LIDAR_NTYPE                       !OPAQ
+    do ip = 1,Npoints                         !OPAQ
+     if(stlidar%cldtype(ip,k).eq.R_UNDEF)then !OPAQ
+       stlidar%cldtype(ip,k)=missing_val      !OPAQ
+     endif                                    !OPAQ
+    enddo                                     !OPAQ
+   enddo                                      !OPAQ
+
 ! AI 11 / 2015
 
@@ -144,4 +165,5 @@
    where(stlidar%cldlayerphase == R_UNDEF) stlidar%cldlayerphase = missing_val
    where(stlidar%lidarcldphase == R_UNDEF) stlidar%lidarcldphase = missing_val
+   where(stlidar%lidarcldtype == R_UNDEF) stlidar%lidarcldtype = missing_val   !OPAQ
    where(stlidar%lidarcldtmp == R_UNDEF) stlidar%lidarcldtmp = missing_val
    
@@ -176,13 +198,28 @@
    CALL histwrite3d_cosp(o_clcalipsotmpun,stlidar%lidarcldtmp(:,:,4),nverttemp)
 
+   CALL histwrite2d_cosp(o_clopaquecalipso,stlidar%cldtype(:,1))               !OPAQ
+   CALL histwrite2d_cosp(o_clthincalipso,stlidar%cldtype(:,2))                 !OPAQ
+   CALL histwrite2d_cosp(o_clzopaquecalipso,stlidar%cldtype(:,3))              !OPAQ
+
+   CALL histwrite3d_cosp(o_clcalipsoopaque,stlidar%lidarcldtype(:,:,1),nvert)  !OPAQ
+   CALL histwrite3d_cosp(o_clcalipsothin,stlidar%lidarcldtype(:,:,2),nvert)    !OPAQ
+   CALL histwrite3d_cosp(o_clcalipsozopaque,stlidar%lidarcldtype(:,:,3),nvert) !OPAQ
+   CALL histwrite3d_cosp(o_clcalipsoopacity,stlidar%lidarcldtype(:,:,4),nvert) !OPAQ
+
+   CALL histwrite3d_cosp(o_proftemp,stlidar%proftemp,nvert)                    !TIBO
+
 #ifdef CPP_XIOS
    CALL histwrite4d_cosp(o_cfad_lidarsr532,stlidar%cfad_sr)
+   CALL histwrite4d_cosp(o_profSR,stlidar%profSR)                              !TIBO
 #else
    do icl=1,SR_BINS
       CALL histwrite3d_cosp(o_cfad_lidarsr532,stlidar%cfad_sr(:,icl,:),nvert,icl)
    enddo
-#endif
-
-   CALL histwrite3d_cosp(o_parasol_refl,stlidar%parasolrefl,nvertp)
+   do icl=1,Ncolumns                                                              !TIBO
+      CALL histwrite3d_cosp(o_profSR,stlidar%profSR(:,icl,:),nvert,icl)           !TIBO
+   enddo                                                                          !TIBO
+     CALL histwrite3d_cosp(o_parasol_refl,stlidar%parasolrefl,nvertp)
+#endif
+
    do k=1,PARASOL_NREFL
     do ip=1, Npoints
Index: /LMDZ5/trunk/libf/phylmd/cosp/cosp_stats.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/cosp_stats.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/cosp_stats.F90	(revision 2955)
@@ -121,4 +121,7 @@
                                            Nlr,vgrid%zl,vgrid%zu,temp_out)
             temp_c(:,:) = temp_out(:,1,:)
+            stlidar%proftemp = temp_c                                     !TIBO
+            where (stlidar%proftemp  < 150.) stlidar%proftemp   = R_UNDEF !TIBO
+            where (stlidar%proftemp  > 350.) stlidar%proftemp   = R_UNDEF !TIBO
 
             call cosp_change_vertical_grid(Npoints,1,Nlevels,gbx%zlev,gbx%zlev_half,ph_in, &
@@ -131,7 +134,8 @@
                             ,LIDAR_UNDEF,ok_lidar_cfad &
                             ,stlidar%cfad_sr,stlidar%srbval &
-                            ,LIDAR_NCAT,stlidar%lidarcld,stlidar%lidarcldphase &
-                            ,stlidar%cldlayer,stlidar%cldlayerphase,stlidar%lidarcldtmp &
-                            ,stlidar%parasolrefl)
+                            ,LIDAR_NCAT,LIDAR_NTYPE,stlidar%lidarcld,stlidar%lidarcldtype & !OPAQ
+                            ,stlidar%lidarcldphase,stlidar%cldlayer,stlidar%cldtype &       !OPAQ
+                            ,stlidar%cldlayerphase,stlidar%lidarcldtmp &                    !OPAQ
+                            ,stlidar%parasolrefl,vgrid%z,stlidar%profSR)                    !OPAQ !TIBO
         endif
 
@@ -140,5 +144,5 @@
                                     temp_c,betatot_out,betaperptot_out,betamol_c,Ze_out, &
                                     stradar%lidar_only_freq_cloud,stradar%radar_lidar_tcc)
-        deallocate(temp_out,temp_c,betaperptot_out)
+        deallocate(temp_in,temp_out,temp_c,betaperptot_out) !TIBO +temp_in
 
         ! Deallocate arrays at coarse resolution
@@ -154,6 +158,8 @@
                         ,LIDAR_UNDEF,ok_lidar_cfad &
                         ,stlidar%cfad_sr,stlidar%srbval &
-                        ,LIDAR_NCAT,stlidar%lidarcld,stlidar%lidarcldphase,stlidar%cldlayer,stlidar%cldlayerphase &
-                        ,stlidar%lidarcldtmp,stlidar%parasolrefl)
+                        ,LIDAR_NCAT,LIDAR_NTYPE,stlidar%lidarcld,stlidar%lidarcldtype & !OPAQ
+                        ,stlidar%lidarcldphase,stlidar%cldlayer,stlidar%cldtype &       !OPAQ
+                        ,stlidar%cldlayerphase,stlidar%lidarcldtmp &                    !OPAQ
+                        ,stlidar%parasolrefl,vgrid%z,stlidar%profSR)                    !OPAQ !TIBO
         !++++++++++++ Lidar-only cloud amount and lidar&radar total cloud mount ++++++++++++++++
         if (cfg%Lradar_sim.and.cfg%Llidar_sim) call cosp_lidar_only_cloud(Npoints,Ncolumns,Nlr, &
@@ -163,6 +169,9 @@
    ! Replace undef
    where (stlidar%cfad_sr   == LIDAR_UNDEF) stlidar%cfad_sr   = R_UNDEF
+   where (stlidar%profSR   == LIDAR_UNDEF) stlidar%profSR   = R_UNDEF !TIBO
    where (stlidar%lidarcld  == LIDAR_UNDEF) stlidar%lidarcld  = R_UNDEF
+   where (stlidar%lidarcldtype  == LIDAR_UNDEF) stlidar%lidarcldtype  = R_UNDEF !OPAQ
    where (stlidar%cldlayer  == LIDAR_UNDEF) stlidar%cldlayer  = R_UNDEF
+   where (stlidar%cldtype  == LIDAR_UNDEF) stlidar%cldtype  = R_UNDEF           !OPAQ
    where (stlidar%parasolrefl == LIDAR_UNDEF) stlidar%parasolrefl = R_UNDEF
    where (stlidar%cldlayerphase  == LIDAR_UNDEF) stlidar%cldlayerphase  = R_UNDEF
Index: /LMDZ5/trunk/libf/phylmd/cosp/cosp_types.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/cosp_types.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/cosp_types.F90	(revision 2955)
@@ -51,5 +51,7 @@
                 Lcltmodis,Lclwmodis,Lclimodis,Lclhmodis,Lclmmodis,Lcllmodis,Ltautmodis,Ltauwmodis,Ltauimodis,Ltautlogmodis, &
                 Ltauwlogmodis,Ltauilogmodis,Lreffclwmodis,Lreffclimodis,Lpctmodis,Llwpmodis, &
-                Liwpmodis,Lclmodis,Lcrimodis,Lcrlmodis
+                Liwpmodis,Lclmodis,Lcrimodis,Lcrlmodis,Lclopaquecalipso,Lclthincalipso, & !OPAQ (2)
+                Lclzopaquecalipso,Lclcalipsoopaque,Lclcalipsothin,Lclcalipsozopaque,    & !OPAQ (4)
+                Lclcalipsoopacity,LprofSR,Lproftemp                                       !OPAQ (1) !TIBO (2)
 
      character(len=32) :: out_list(N_OUT_LIST)
@@ -199,9 +201,14 @@
     ! Arrays with dimensions (Npoints,Nlevels)
     real, dimension(:,:),pointer :: lidarcld    ! 3D "lidar" cloud fraction 
+    real, dimension(:,:),pointer :: proftemp    ! Temperature profiles 40 levs !TIBO 
     ! Arrays with dimensions (Npoints,LIDAR_NCAT)
     real, dimension(:,:),pointer :: cldlayer      ! low, mid, high-level, total lidar cloud cover
-   ! Arrays with dimensions (Npoints,Nlevels,Nphase)
+    ! Arrays with dimensions (Npoints,LIDAR_NTYPE)                                               !OPAQ
+    real, dimension(:,:),pointer :: cldtype       ! opaque and thin cloud covers, z_opaque       !OPAQ
+    ! Arrays with dimensions (Npoints,Nlevels,Nphase)
     real, dimension(:,:,:),pointer :: lidarcldphase    ! 3D "lidar" phase cloud fraction 
-     ! Arrays with dimensions (Npoints,LIDAR_NCAT,Nphase)
+    ! Arrays with dimensions (Npoints,Nlevels,LIDAR_NTYPE+1)                                     !OPAQ
+    real, dimension(:,:,:),pointer :: lidarcldtype     ! 3D "lidar" OPAQ type fraction + opacity !OPAQ 
+    ! Arrays with dimensions (Npoints,LIDAR_NCAT,Nphase)
     real, dimension(:,:,:),pointer :: cldlayerphase      ! low, mid, high-level lidar phase cloud cover
     ! Arrays with dimensions (Npoints,Ntemps,Nphase)
@@ -209,4 +216,8 @@
     ! Arrays with dimensions (Npoints,PARASOL_NREFL)
     real, dimension(:,:),pointer :: parasolrefl   ! mean parasol reflectance
+!    ! Arrays with dimensions (Npoints,Ncolumns,Nlevels)                     !TIBO
+!    real, dimension(:,:,:),pointer :: profSR      ! subcolumns for each day !TIBO 
+    ! Arrays with dimensions (Npoints,Nlevels,Ncolumns)                     !TIBO2
+    real, dimension(:,:,:),pointer :: profSR      ! subcolumns for each day !TIBO2 
 
   END TYPE COSP_LIDARSTATS
@@ -803,4 +814,7 @@
     allocate(x%lidarcldphase(i,k,6),x%lidarcldtmp(i,LIDAR_NTEMP,5),&
              x%cldlayerphase(i,LIDAR_NCAT,6))
+    allocate(x%lidarcldtype(i,k,LIDAR_NTYPE+1),x%cldtype(i,LIDAR_NTYPE)) !OPAQ
+!    allocate(x%profSR(i,j,k),x%proftemp(i,k))                            !TIBO
+    allocate(x%profSR(i,k,j),x%proftemp(i,k))                            !TIBO2
     ! --- Initialise to zero ---
     x%srbval    = 0.0
@@ -812,4 +826,8 @@
     x%cldlayerphase  = 0.0
     x%lidarcldtmp  = 0.0
+    x%lidarcldtype  = 0.0 !OPAQ
+    x%cldtype  = 0.0      !OPAQ
+    x%profSR   = 0.0      !TIBO
+    x%proftemp = 0.0      !TIBO
 
    END SUBROUTINE CONSTRUCT_COSP_LIDARSTATS
@@ -823,4 +841,6 @@
     deallocate(x%srbval, x%cfad_sr, x%lidarcld, x%cldlayer, x%parasolrefl)
     deallocate(x%cldlayerphase, x%lidarcldtmp, x%lidarcldphase)
+    deallocate(x%lidarcldtype, x%cldtype) !OPAQ
+    deallocate(x%profSR, x%proftemp)      !TIBO
   END SUBROUTINE FREE_COSP_LIDARSTATS
 
@@ -1367,4 +1387,8 @@
     y%cldlayerphase(iy(1):iy(2),:,:)  = x%cldlayerphase(ix(1):ix(2),:,:)
     y%lidarcldtmp(iy(1):iy(2),:,:)    = x%lidarcldtmp(ix(1):ix(2),:,:)
+    y%lidarcldtype(iy(1):iy(2),:,:)    = x%lidarcldtype(ix(1):ix(2),:,:) !OPAQ
+    y%cldtype(iy(1):iy(2),:)           = x%cldtype(ix(1):ix(2),:)        !OPAQ
+    y%profSR(iy(1):iy(2),:,:)          = x%profSR(ix(1):ix(2),:,:)       !TIBO
+    y%proftemp(iy(1):iy(2),:)          = x%proftemp(ix(1):ix(2),:)       !TIBO
 END SUBROUTINE COSP_LIDARSTATS_CPSECTION
 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -1599,12 +1623,21 @@
     ! Arrays with dimensions (Npoints,SR_BINS,Nlevels)
     print *, 'shape(x%cfad_sr): ',shape(x%cfad_sr)
+!    ! Arrays with dimensions (Npoints,Ncolumns,Nlevels) !TIBO
+!    print *, 'shape(x%profSR): ',shape(x%profSR)        !TIBO
+    ! Arrays with dimensions (Npoints,Nlevels,Ncolumns) !TIBO2
+    print *, 'shape(x%profSR): ',shape(x%profSR)        !TIBO2
     ! Arrays with dimensions (Npoints,Nlevels)
     print *, 'shape(x%lidarcld): ',shape(x%lidarcld)
+    print *, 'shape(x%proftemp): ',shape(x%proftemp)    !TIBO
     ! Arrays with dimensions (Npoints,LIDAR_NCAT)
     print *, 'shape(x%cldlayer): ',shape(x%cldlayer)
+    ! Arrays with dimensions (Npoints,LIDAR_NTYPE)            !OPAQ
+    print *, 'shape(x%cldtype): ',shape(x%cldtype)            !OPAQ
     ! Arrays with dimensions (Npoints,PARASOL_NREFL)
     print *, 'shape(x%parasolrefl): ',shape(x%parasolrefl)
      ! Arrays with dimensions (Npoints,Nlevels,Nphase)
     print *, 'shape(x%lidarcldphase): ',shape(x%lidarcldphase)
+     ! Arrays with dimensions (Npoints,Nlevels,LIDAR_NTYPE+1) !OPAQ
+    print *, 'shape(x%lidarcldtype): ',shape(x%lidarcldtype)  !OPAQ
      ! Arrays with dimensions (Npoints,LIDAR_NCAT,Nphase)
     print *, 'shape(x%cldlayerphase): ',shape(x%cldlayerphase)
Index: /LMDZ5/trunk/libf/phylmd/cosp/lmd_ipsl_stats.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/lmd_ipsl_stats.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/lmd_ipsl_stats.F90	(revision 2955)
@@ -36,6 +36,6 @@
       SUBROUTINE diag_lidar(npoints,ncol,llm,max_bin,nrefl &
                   ,tmp,pnorm,pnorm_perp,pmol,refl,land,pplay,undef,ok_lidar_cfad &
-                  ,cfad2,srbval,ncat,lidarcld,lidarcldphase,cldlayer,cldlayerphase &
-                  ,lidarcldtmp,parasolrefl)
+                  ,cfad2,srbval,ncat,ntype,lidarcld,lidarcldtype,lidarcldphase,cldlayer & !OPAQ
+                  ,cldtype,cldlayerphase,lidarcldtmp,parasolrefl,vgrid_z,profSR)          !OPAQ !TIBO
 !
 ! -----------------------------------------------------------------------------------
@@ -82,4 +82,5 @@
       integer max_bin               ! nb of bins for SR CFADs
       integer ncat                  ! nb of cloud layer types (low,mid,high,total)
+      integer ntype                 ! nb of OPAQ products (opaque and thin clouds, z_opaque) !OPAQ
       integer nrefl                 ! nb of solar zenith angles for parasol reflectances
 
@@ -93,13 +94,18 @@
       real tmp(npoints,llm)         ! temp at each levels
       real pnorm_perp(npoints,ncol,llm)  ! lidar perpendicular ATB
+      real vgrid_z(llm)             ! mid-level altitude of the output vertical grid         !OPAQ
 
 ! c outputs :
       real lidarcld(npoints,llm)     ! 3D "lidar" cloud fraction
+      real lidarcldtype(npoints,llm,ntype+1)   ! 3D "lidar" OPAQ type fraction + opacity     !OPAQ
       real sub(npoints,llm)     ! 3D "lidar" indice
       real cldlayer(npoints,ncat)    ! "lidar" cloud layer fraction (low, mid, high, total)
+      real cldtype(npoints,ntype)  ! "lidar" OPAQ type covers (opaque/thin cloud + z_opaque) !OPAQ
 
       real cfad2(npoints,max_bin,llm) ! CFADs of SR
       real srbval(max_bin)           ! SR bins in CFADs
       real parasolrefl(npoints,nrefl)! grid-averaged parasol reflectance
+!     real profSR(npoints,ncol,llm)  ! tableau avec les subcolumns SR !TIBO
+      real profSR(npoints,llm,ncol)  ! tableau avec les subcolumns SR !TIBO2
 
 ! c threshold for cloud detection :
@@ -110,4 +116,5 @@
       real S_att
       parameter (S_att = 0.01)
+!      parameter (S_att = 0.06)  !OPAQ ! Threshold for "surface detection" equivalent
 
 ! c local variables :
@@ -147,4 +154,6 @@
         end where
          x3d(:,ic,:) = x3d_c
+!	profSR(:,ic,:) = x3d(:,ic,:) !TIBO
+	profSR(:,:,ic) = x3d(:,ic,:) !TIBO2
       enddo
 
@@ -157,4 +166,7 @@
               tmp,x3d,pnorm,pnorm_perp,pplay, S_att,S_cld,S_cld_att,undef,lidarcld, &
               cldlayer,lidarcldphase,sub,cldlayerphase,lidarcldtmp)
+
+    CALL COSP_OPAQ(npoints,ncol,llm,ntype,x3d,S_cld,undef,lidarcldtype,            & !OPAQ
+                   cldtype,vgrid_z)                                                  !OPAQ
 
 ! c -------------------------------------------------------
@@ -983,4 +995,187 @@
 ! ---------------------------------------------------------------
 
+! BEGINNING OF OPAQ CHANGES
+    ! ####################################################################################
+    ! SUBROUTINE cosp_opaq
+    ! Conventions: Ntype must be equal to 3 (opaque cloud, thin cloud, z_opaque)
+    ! ####################################################################################
+    SUBROUTINE COSP_OPAQ(Npoints,Ncolumns,Nlevels,Ntype,x,S_cld,undef,lidarcldtype,   &
+                         cldtype,vgrid_z)
+
+      IMPLICIT NONE
+! Input arguments
+      integer Npoints,Ncolumns,Nlevels,Ntype
+      real x(Npoints,Ncolumns,Nlevels)
+      real S_cld
+      real undef
+      real vgrid_z(Nlevels)
+! Output :
+      real lidarcldtype(Npoints,Nlevels,Ntype+1) ! 3D "lidar" OPAQ type + opacity fraction
+      real cldtype(Npoints,Ntype)              ! opaque and thin cloud covers, z_opaque
+! Local variables
+      integer ip, k, iz, ic, zopac
+      real p1
+      real cldy(Npoints,Ncolumns,Nlevels)
+      real cldyopaq(Npoints,Ncolumns,Nlevels)
+      real srok(Npoints,Ncolumns,Nlevels)
+      real srokopaq(Npoints,Ncolumns,Nlevels)
+      real cldlay(Npoints,Ncolumns,Ntype+1)  ! opaque, thin, z_opaque and all cloud cover
+      real nsublay(Npoints,Ncolumns,Ntype+1) ! opaque, thin, z_opaque and all cloud cover
+      real nsublayer(Npoints,Ntype)
+      real nsub(Npoints,Nlevels)
+      real nsubopaq(Npoints,Nlevels)
+      real S_att_opaq
+      real S_att
+  
+    ! ####################################################################################
+	! 1) Initialize    
+    ! ####################################################################################
+    cldtype               = 0.0
+    lidarcldtype          = 0.0
+    nsub                  = 0.0
+    nsubopaq              = 0.0
+    cldlay                = 0.0
+    nsublay               = 0.0
+    nsublayer             = 0.0
+    S_att_opaq            = 0.06 ! Fully Attenuated threshold, from Guzman et al. 2017, JGR-A
+    S_att                 = 0.01 
+
+    ! ####################################################################################
+    ! 2) Cloud detection and Fully attenuated layer detection
+    ! ####################################################################################
+    do k=1,Nlevels
+       ! Cloud detection at subgrid-scale:
+       where ( (x(:,:,k) .gt. S_cld) .and. (x(:,:,k) .ne. undef) )
+          cldy(:,:,k)=1.0
+       elsewhere
+          cldy(:,:,k)=0.0
+       endwhere
+       ! Fully attenuated layer detection at subgrid-scale:
+       where ( (x(:,:,k) .gt. 0.0) .and. (x(:,:,k) .lt. S_att_opaq) .and. (x(:,:,k) .ne. undef) )
+          cldyopaq(:,:,k)=1.0
+       elsewhere
+          cldyopaq(:,:,k)=0.0
+       endwhere
+
+       ! Number of useful sub-column layers:
+       where ( (x(:,:,k) .gt. S_att) .and. (x(:,:,k) .ne. undef) )
+          srok(:,:,k)=1.0
+       elsewhere
+          srok(:,:,k)=0.0
+       endwhere
+       ! Number of useful sub-columns layers for z_opaque 3D fraction:
+       where ( (x(:,:,k) .gt. 0.0) .and. (x(:,:,k) .ne. undef) )
+          srokopaq(:,:,k)=1.0
+       elsewhere
+          srokopaq(:,:,k)=0.0
+       endwhere
+    enddo
+
+    ! ####################################################################################
+    ! 3) Grid-box 3D OPAQ product fraction and cloud type cover (opaque/thin) + mean z_opaque
+    ! ####################################################################################
+
+    do k= Nlevels,1,-1
+       do ic = 1, Ncolumns
+          do ip = 1, Npoints
+
+             cldlay(ip,ic,1)   = MAX(cldlay(ip,ic,1),cldyopaq(ip,ic,k)) ! Opaque clouds
+             cldlay(ip,ic,4)   = MAX(cldlay(ip,ic,4),cldy(ip,ic,k))     ! All clouds
+
+             nsublay(ip,ic,1)  = MAX(nsublay(ip,ic,1),srok(ip,ic,k))
+             nsublay(ip,ic,2)  = MAX(nsublay(ip,ic,2),srok(ip,ic,k))
+!             nsublay(ip,ic,4)  = MAX(nsublay(ip,ic,4),srok(ip,ic,k))
+             nsub(ip,k)        = nsub(ip,k) + srok(ip,ic,k)
+             nsubopaq(ip,k)    = nsubopaq(ip,k) + srokopaq(ip,ic,k)
+
+          enddo
+       enddo
+    enddo   
+
+! OPAQ variables
+     do ic = 1, Ncolumns
+        do ip = 1, Npoints
+
+     ! Declaring non-opaque cloudy profiles as thin cloud profiles
+	   if ( (cldlay(ip,ic,4) .eq. 1.0) .and. (cldlay(ip,ic,1) .eq. 0.0) ) then
+	      cldlay(ip,ic,2)  =  1.0
+ 	   endif
+
+     ! Filling in 3D and 2D variables
+
+     ! Opaque cloud profiles
+	   if ( cldlay(ip,ic,1) .eq. 1.0 ) then
+	      zopac = 0.0
+	      do k=2,Nlevels
+     ! Declaring opaque cloud fraction and z_opaque altitude for 3D and 2D variables
+	         if ( (cldy(ip,ic,k) .eq. 1.0) .and. (zopac .eq. 0.0) ) then
+		    lidarcldtype(ip,k-1,3) = lidarcldtype(ip,k-1,3) + 1.0
+		    cldlay(ip,ic,3)        = vgrid_z(k-1) !z_opaque altitude
+		    nsublay(ip,ic,3)       = 1.0
+		    zopac = 1.0
+		 endif
+	         if ( cldy(ip,ic,k) .eq. 1.0 ) then
+		    lidarcldtype(ip,k,1)   = lidarcldtype(ip,k,1) + 1.0
+                 endif
+	      enddo
+	   endif
+
+     ! Thin cloud profiles
+	   if ( cldlay(ip,ic,2) .eq. 1.0 ) then
+	      do k=1,Nlevels
+     ! Declaring thin cloud fraction for 3D variable
+                 if ( cldy(ip,ic,k) .eq. 1.0 ) then
+                    lidarcldtype(ip,k,2) = lidarcldtype(ip,k,2) + 1.0
+                 endif
+	      enddo
+           endif
+
+       enddo
+    enddo   
+
+    ! 3D cloud types fraction (opaque=1 and thin=2)
+    where ( nsub(:,:) .gt. 0.0 )
+       lidarcldtype(:,:,1) = lidarcldtype(:,:,1)/nsub(:,:)
+       lidarcldtype(:,:,2) = lidarcldtype(:,:,2)/nsub(:,:)
+    elsewhere
+       lidarcldtype(:,:,1) = undef
+       lidarcldtype(:,:,2) = undef
+    endwhere
+    ! 3D z_opaque fraction (=3)
+    where ( nsubopaq(:,:) .gt. 0.0 )
+       lidarcldtype(:,:,3) = lidarcldtype(:,:,3)/nsubopaq(:,:)
+    elsewhere
+       lidarcldtype(:,:,3) = undef
+    endwhere
+    ! 3D opacity fraction (=4) !Summing z_opaque fraction from TOA(k=Nlevels) to SFC(k=1)
+       lidarcldtype(:,Nlevels,4) = lidarcldtype(:,Nlevels,3)
+    do ip = 1, Npoints
+     	do k = Nlevels-1, 1, -1
+           if ( lidarcldtype(ip,k,3) .ne. undef ) then
+	      lidarcldtype(ip,k,4) = lidarcldtype(ip,k+1,4) + lidarcldtype(ip,k,3)
+           endif
+	enddo
+    enddo
+    where ( nsubopaq(:,:) .eq. 0.0 )
+       lidarcldtype(:,:,4) = undef
+    endwhere
+
+    ! Layered cloud types (opaque, thin and z_opaque 2D variables)
+
+    do iz = 1, Ntype
+       do ic = 1, Ncolumns
+          cldtype(:,iz)   = cldtype(:,iz)   + cldlay(:,ic,iz)
+          nsublayer(:,iz) = nsublayer(:,iz) + nsublay(:,ic,iz)
+       enddo
+    enddo
+    where (nsublayer(:,:) .gt. 0.0)
+       cldtype(:,:) = cldtype(:,:)/nsublayer(:,:)
+    elsewhere
+       cldtype(:,:) = undef
+    endwhere
+
+  END SUBROUTINE COSP_OPAQ
+! END OF OPAQ CHANGES
+
 
 END MODULE MOD_LMD_IPSL_STATS
Index: /LMDZ5/trunk/libf/phylmd/cosp/read_cosp_output_nl.F90
===================================================================
--- /LMDZ5/trunk/libf/phylmd/cosp/read_cosp_output_nl.F90	(revision 2954)
+++ /LMDZ5/trunk/libf/phylmd/cosp/read_cosp_output_nl.F90	(revision 2955)
@@ -41,5 +41,7 @@
              Lcltmodis,Lclwmodis,Lclimodis,Lclhmodis,Lclmmodis,Lcllmodis,Ltautmodis,Ltauwmodis,Ltauimodis,Ltautlogmodis, &
              Ltauwlogmodis,Ltauilogmodis,Lreffclwmodis,Lreffclimodis,Lpctmodis,Llwpmodis, &
-             Liwpmodis,Lclmodis,Lcrimodis,Lcrlmodis
+             Liwpmodis,Lclmodis,Lcrimodis,Lcrlmodis,Lclopaquecalipso,Lclthincalipso,      &           !OPAQ (2)
+             Lclzopaquecalipso,Lclcalipsoopaque,Lclcalipsothin,Lclcalipsozopaque,Lclcalipsoopacity, & !OPAQ (5)
+             LprofSR,Lproftemp                                                                        !TIBO (2)
 
   namelist/COSP_OUTPUT/Lradar_sim,Llidar_sim,Lisccp_sim,Lmodis_sim,Lmisr_sim,Lrttov_sim, &
@@ -57,5 +59,7 @@
              Lcltmodis,Lclwmodis,Lclimodis,Lclhmodis,Lclmmodis,Lcllmodis,Ltautmodis,Ltauwmodis,Ltauimodis,Ltautlogmodis, &
              Ltauwlogmodis,Ltauilogmodis,Lreffclwmodis,Lreffclimodis,Lpctmodis,Llwpmodis, &
-             Liwpmodis,Lclmodis,Lcrimodis,Lcrlmodis
+             Liwpmodis,Lclmodis,Lcrimodis,Lcrlmodis,Lclopaquecalipso,Lclthincalipso,      &           !OPAQ (2)
+             Lclzopaquecalipso,Lclcalipsoopaque,Lclcalipsothin,Lclcalipsozopaque,Lclcalipsoopacity, & !OPAQ (5)
+             LprofSR,Lproftemp                                                                        !TIBO (2)
    
   do i=1,N_OUT_LIST
@@ -139,4 +143,13 @@
   CALL bcast(Lcrimodis)
   CALL bcast(Lcrlmodis)
+  CALL bcast(Lclopaquecalipso)  !OPAQ
+  CALL bcast(Lclthincalipso)    !OPAQ
+  CALL bcast(Lclzopaquecalipso) !OPAQ
+  CALL bcast(Lclcalipsoopaque)  !OPAQ
+  CALL bcast(Lclcalipsothin)    !OPAQ
+  CALL bcast(Lclcalipsozopaque) !OPAQ
+  CALL bcast(Lclcalipsoopacity) !OPAQ
+  CALL bcast(LprofSR)           !TIBO
+  CALL bcast(Lproftemp)         !TIBO
 
 !$OMP BARRIER
@@ -185,5 +198,15 @@
     Lclmcalipsoun      = .false.
     Lcltcalipsoun      = .false.
-  endif
+    Lclopaquecalipso   = .false. !OPAQ
+    Lclthincalipso     = .false. !OPAQ
+    Lclzopaquecalipso  = .false. !OPAQ
+    Lclcalipsoopaque   = .false. !OPAQ
+    Lclcalipsothin     = .false. !OPAQ
+    Lclcalipsozopaque  = .false. !OPAQ
+    Lclcalipsoopacity  = .false. !OPAQ
+    LprofSR            = .false. !TIBO
+    Lproftemp          = .false. !TIBO
+  endif
+
   if (.not.Lisccp_sim) then
     Lalbisccp       = .false.
@@ -390,4 +413,23 @@
   i = i+1
   if (Lcrlmodis)         cfg%out_list(i) = 'crlmodis'
+
+  i = i+1                                                            !OPAQ
+  if (Lclopaquecalipso)         cfg%out_list(i) = 'clopaquecalipso'  !OPAQ
+  i = i+1                                                            !OPAQ
+  if (Lclthincalipso)           cfg%out_list(i) = 'clthincalipso'    !OPAQ
+  i = i+1                                                            !OPAQ
+  if (Lclzopaquecalipso)        cfg%out_list(i) = 'clzopaquecalipso' !OPAQ
+  i = i+1                                                            !OPAQ
+  if (Lclcalipsoopaque)         cfg%out_list(i) = 'clcalipsoopaque'  !OPAQ
+  i = i+1                                                            !OPAQ
+  if (Lclcalipsothin)           cfg%out_list(i) = 'clcalipsothin'    !OPAQ
+  i = i+1                                                            !OPAQ
+  if (Lclcalipsozopaque)        cfg%out_list(i) = 'clcalipsozopaque' !OPAQ
+  i = i+1                                                            !OPAQ
+  if (Lclcalipsoopacity)        cfg%out_list(i) = 'clcalipsoopacity' !OPAQ
+  i = i+1                                                            !TIBO
+  if (LprofSR)                  cfg%out_list(i) = 'profSR'           !TIBO
+  i = i+1                                                            !TIBO
+  if (Lproftemp)                cfg%out_list(i) = 'proftemp'         !TIBO
     
   if (i /= N_OUT_LIST) then
@@ -441,4 +483,13 @@
   cfg%Lcltlidarradar = Lcltlidarradar
   cfg%LparasolRefl = LparasolRefl
+  cfg%Lclopaquecalipso  = Lclopaquecalipso  !OPAQ
+  cfg%Lclthincalipso    = Lclthincalipso    !OPAQ
+  cfg%Lclzopaquecalipso = Lclzopaquecalipso !OPAQ
+  cfg%Lclcalipsoopaque  = Lclcalipsoopaque  !OPAQ
+  cfg%Lclcalipsothin    = Lclcalipsothin    !OPAQ
+  cfg%Lclcalipsozopaque = Lclcalipsozopaque !OPAQ
+  cfg%Lclcalipsoopacity = Lclcalipsoopacity !OPAQ
+  cfg%LprofSR           = LprofSR           !TIBO
+  cfg%Lproftemp         = Lproftemp         !TIBO
   ! MISR simulator  
   cfg%LclMISR = LclMISR
