<p><b>ringler@lanl.gov</b> 2011-10-14 11:39:11 -0600 (Fri, 14 Oct 2011)</p><p><br>
module_cullLoops is working in beta mode.<br>
it successfully removed all inland seas in the 10242 mesh<br>
the default in basin.F is true, i.e. inland seas are removed<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/basin/src/basin.F
===================================================================
--- branches/ocean_projects/basin/src/basin.F        2011-10-14 15:44:37 UTC (rev 1076)
+++ branches/ocean_projects/basin/src/basin.F        2011-10-14 17:39:11 UTC (rev 1077)
@@ -5,7 +5,7 @@
use read_TS
use write_netcdf
use utilities
-!use cullLoops
+use cullLoops
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Program: basin.F
@@ -83,6 +83,7 @@
!real*8, parameter :: sphere_radius = 0.0
logical, parameter :: real_bathymetry=.true.
+logical, parameter :: eliminate_inland_seas=.true.
logical, parameter :: l_woce = .true.
@@ -1193,10 +1194,12 @@
if(flag) kmt(iCell)=0
enddo
-!call eliminateLoops(nCells,nEdges,nVertices,maxEdges,vertexDegree, &
- !nEdgesOnCell, cellsOnCell, verticesOnEdge, cellsOnVertex, edgesOnCell, lonCell, latCell, &
- !xCell, yCell, zCell, xEdge, yEdge, zEdge, xVertex, yVertex, zVertex, &
- !KMT)
+if(eliminate_inland_seas) then
+call eliminateLoops(nCells,nEdges,nVertices,maxEdges,vertexDegree, &
+ nEdgesOnCell, cellsOnCell, verticesOnEdge, cellsOnVertex, edgesOnCell, lonCell, latCell, &
+ xCell, yCell, zCell, xEdge, yEdge, zEdge, xVertex, yVertex, zVertex, &
+ KMT)
+endif
if(real_bathymetry) then
where(kmt.eq.1) kmt=3
Modified: branches/ocean_projects/basin/src/module_cullLoops.F
===================================================================
--- branches/ocean_projects/basin/src/module_cullLoops.F        2011-10-14 15:44:37 UTC (rev 1076)
+++ branches/ocean_projects/basin/src/module_cullLoops.F        2011-10-14 17:39:11 UTC (rev 1077)
@@ -25,14 +25,17 @@
integer, intent(inout) :: KMT(ncells)
! local workspace
- integer :: iCell, jCell, oCell, lCell, iEdge, i, kCell, iSharedEdge, iStartEdge, iSave
+ integer :: iCell, jCell, oCell, lCell, iEdge, i, kCell, iSharedEdge, iStartEdge, iSave, iSweep
integer :: iEdgeCounter, nEdgesInLoop(nCells), iCellAhead, LeftTurns, RightTurns
- logical :: connected, atBoundary
- real :: lat, rlat, rCenter(3), s(3), t(3), q(3), rCross
+ logical :: connected, atBoundary, moveSouth, moveEast, atGrenwich
+ real :: lat, rlat, rlon, rCenter(3), s(3), t(3), q(3), rCross, mylon, mylat, pi
+ pi = 4.0*atan(1.0)
+
! we loop over all cells and count the number of edges in the loop containing iCell
! there is no coupling between iCell, so this can be inside an openMP directive
iCellMask(:) = 0
+ moveSouth = .true.
do iCell=1,nCells
write(6,*) 'working on : ',iCell, KMT(iCell)
@@ -66,7 +69,9 @@
! choose the next cell to be the one with the most southern latitude
! this jCell will only be used if atBoundary=.false., thus jCell must be an ocean cell
+ if(moveSouth) then
rlat = 10000.0
+ mylat = latCell(jCell)
do i=1,nEdgesOnCell(jCell)
kCell = cellsOnCell(i,jCell)
if(latCell(kCell).lt.rlat) then
@@ -75,8 +80,9 @@
endif
enddo
jCell = iSave
+ endif
- if(.not.atBoundary) write(6,*) ' pushing on to the south : ', jCell
+ if(moveSouth.and..not.atBoundary) write(6,*) ' pushing on to the south : ', jCell
enddo ! .not.atBoundary
@@ -151,10 +157,9 @@
write(6,*)
write(6,*) ' edges and cross ', iEdgeCounter, rCross, LeftTurns, RightTurns
- if(rCross.gt.0.0) write(6,*) ' inland sea was found ', LeftTurns-RightTurns
write(6,*)
- if(rCross.gt.0.0) then
+ if(LeftTurns-RightTurns.gt.0.and.rCross.gt.0.0) then
iCellMask(iCell) = 1
write(50,11) iCell, lonCell(iCell), latCell(iCell)
11 format(i5,2f12.4)
@@ -163,6 +168,18 @@
enddo ! iCell
! cull all inland seas
+ do iSweep=1,nCells/10
+ write(6,*) iSweep
+ do iCell=1,nCells
+ if(iCellMask(iCell).eq.1) then
+ do i=1,nEdgesOnCell(iCell)
+ kCell=cellsOnCell(i,iCell)
+ if(KMT(kCell).gt.0) iCellMask(kCell)=1
+ enddo
+ endif
+ enddo
+ enddo
+
write(6,*) ' total cells culled ', sum(iCellMask)
where(iCellMask(:).eq.1) KMT(:)=0
</font>
</pre>