<p><b>mhoffman@lanl.gov</b> 2013-02-04 17:40:53 -0700 (Mon, 04 Feb 2013)</p><p>BRANCH COMMIT - land ice<br>
<br>
Modifying the area diagnostic so it only reports cells with ice!<br>
</p><hr noshade><pre><font color="gray">Modified: branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_global_diagnostics.F
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_global_diagnostics.F        2013-02-04 22:46:17 UTC (rev 2426)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_global_diagnostics.F        2013-02-05 00:40:53 UTC (rev 2427)
@@ -1,3 +1,6 @@
+
+#include &quot;mpas_land_ice_mask.inc&quot;
+
 module land_ice_global_diagnostics
 
    use mpas_grid_types
@@ -55,20 +58,21 @@
       real (kind=RKIND), dimension(:,:), pointer :: layerThickness, normalVelocity, &amp;
 !!!                                                    pv_edge, pv_vertex, pv_cell,  &amp;
                                                     layerThicknessVertex, weightsOnEdge
-
+      integer, dimension(:), pointer :: cellMask
       real (kind=RKIND), dimension(:,:,:), pointer :: tracers
 
       real (kind=RKIND), dimension(:), allocatable :: volumeWeightedPotentialEnergyReservoir, averageThickness
       real (kind=RKIND), dimension(:), allocatable :: areaEdge, potentialEnstrophyReservior
 !!!      real (kind=RKIND), dimension(:), allocatable :: h_s_edge
 
+      real (kind=RKIND), dimension(:), allocatable :: iceArea
       real (kind=RKIND), dimension(:,:), allocatable :: cellVolume, cellArea, volumeWeightedPotentialVorticity
       real (kind=RKIND), dimension(:,:), allocatable :: volumeWeightedPotentialEnstrophy, vertexVolume, volumeWeightedKineticEnergy 
       real (kind=RKIND), dimension(:,:), allocatable :: volumeWeightedPotentialEnergy, volumeWeightedPotentialEnergyTopography 
       real (kind=RKIND), dimension(:,:), allocatable :: keTend_CoriolisForce, keTend_PressureGradient 
       real (kind=RKIND), dimension(:,:), allocatable ::peTend_DivThickness, refAreaWeightedSurfaceHeight, refAreaWeightedSurfaceHeight_edge
 
-      real (kind=RKIND) :: sumCellVolume, sumCellArea, sumVertexVolume, sumrefAreaWeightedSurfaceHeight
+      real (kind=RKIND) :: sumCellVolume, sumIceAreaLocal, sumIceArea, sumCellArea, sumVertexVolume, sumrefAreaWeightedSurfaceHeight
 
       real (kind=RKIND) :: globalIceThickness, globalIceMass
       !real (kind=RKIND) :: globalFluidThickness, globalPotentialVorticity, globalPotentialEnstrophy, globalEnergy 
@@ -108,6 +112,7 @@
       areaEdge = dcEdge(1:nEdgesSolve)*dvEdge(1:nEdgesSolve)
       weightsOnEdge     =&gt; grid % weightsOnEdge % array
 
+      cellMask =&gt; state % cellMask % array
       layerThickness =&gt; state % layerThickness % array
       normalVelocity =&gt; state % normalVelocity % array
 !!!      v =&gt; state % v % array
@@ -121,6 +126,7 @@
       ! Step 2
       ! 2. Allocate the array with the correct dimensions.
       allocate(cellVolume(nVertLevels,nCellsSolve))
+      allocate(iceArea(nCells+1))
       allocate(cellArea(nVertLevels,nCellsSolve))
       allocate(refAreaWeightedSurfaceHeight(nVertLevels,nCellsSolve))
       allocate(refAreaWeightedSurfaceHeight_edge(nVertLevels,nEdgesSolve))
@@ -142,6 +148,7 @@
 
 
       cellVolume = 0.0_RKIND
+      iceArea = 0.0_RKIND
       refAreaWeightedSurfaceHeight = 0.0_RKIND
       refAreaWeightedSurfaceHeight_edge = 0.0_RKIND
       vertexVolume = 0.0_RKIND
@@ -162,11 +169,17 @@
       ! Step 3
       ! 3. Fill the array with the data to be integrated.
       !     eg. GlobalFluidThickness = Sum(h dA)/Sum(dA), See below for array filling
+      where ( MASK_IS_ICE( cellMask ) )
+        iceArea = areaCell
+      else where
+        iceArea = 0.0_RKIND
+      end where
+
       do iLevel = 1,nVertLevels
         ! eg. GlobalFluidThickness top (Sum( h dA)) = Sum(cellVolume)
         cellVolume(iLevel,:) = layerThickness(iLevel,1:nCellsSolve)*areaCell(1:nCellsSolve)
         ! eg. GlobalFluidThickness bot (Sum(dA)) = Sum(cellArea)
-        cellArea(iLevel,:) = areaCell(1:nCellsSolve)
+        !cellArea(iLevel,:) = areaCell(1:nCellsSolve)
 !!!        volumeWeightedPotentialVorticity(iLevel,:) = pv_vertex(iLevel,1:nVerticesSolve) &amp;
 !!!                *layerThicknessVertex(iLevel,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
 !!!        volumeWeightedPotentialEnstrophy(iLevel,:) = pv_vertex(iLevel,1:nVerticesSolve) &amp; 
@@ -219,8 +232,11 @@
 
       ! Computing top and bottom of global mass integral
       call land_ice_compute_global_sum(dminfo, nVertLevels, nCellsSolve, cellVolume, sumCellVolume)
-      call land_ice_compute_global_sum(dminfo, nVertLevels, nCellsSolve, cellArea, sumCellArea)
+      sumIceAreaLocal = sum(iceArea(1:nCellsSolve))
+      call mpas_dmpar_sum_real(dminfo, sumIceAreaLocal, sumIceArea)  ! Compute global sum
 
+!      call land_ice_compute_global_sum(dminfo, nVertLevels, nCellsSolve, cellArea, sumCellArea)
+
 !!!      globalKineticEnergyTendency = globalKineticEnergyTendency / sumCellVolume
 !!!      globalPotentialEnergyTendency = globalPotentialEnergyTendency / sumCellVolume
 
@@ -288,7 +304,7 @@
 
          ! Get the precision being used so we don't write more digits than are meaningful
          write(realPrecisionStr,*) PRECISION(sumCellVolume)-1
-         write(fileID,'(1i0, 100es24.'// ADJUSTL(realPrecisionStr) // ')') timeIndex, timeIndex*dt, sumCellArea, sumCellVolume
+         write(fileID,'(1i0, 100es24.'// ADJUSTL(realPrecisionStr) // ')') timeIndex, timeIndex*dt, sumIceArea, sumCellVolume
          close(fileID)
       end if
 

</font>
</pre>