<p><b>ringler@lanl.gov</b> 2011-10-24 15:13:17 -0600 (Mon, 24 Oct 2011)</p><p><br>
1) propagate meshDensity from grid.nc to ocean.nc<br>
(if meshDensity is not present in grid.nc, set to 1.0 in ocean.nc)<br>
<br>
2) add more data to dx output directory<br>
<br>
3) fix the error in shift of 180 degrees<br>
<br>
4) re-implement the raising of KMT based on surrounding cell depths<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-24 20:25:34 UTC (rev 1123)
+++ branches/ocean_projects/basin/src/basin.F 2011-10-24 21:13:17 UTC (rev 1124)
@@ -37,7 +37,7 @@
integer :: time, nCells, nEdges, nVertices
integer :: maxEdges, maxEdges2, TWO, vertexDegree, nVertLevels
integer, allocatable, dimension(:) :: indexToCellID, indexToEdgeID, indexToVertexID
-real, allocatable, dimension(:) :: xCell, yCell, zCell, latCell, lonCell
+real, allocatable, dimension(:) :: xCell, yCell, zCell, latCell, lonCell, meshDensity
real, allocatable, dimension(:) :: xEdge, yEdge, zEdge, latEdge, lonEdge
real, allocatable, dimension(:) :: xVertex, yVertex, zVertex, latVertex, lonVertex
integer, allocatable, dimension(:) :: nEdgesOnCell, nEdgesOnEdge
@@ -105,7 +105,7 @@
integer :: nCellsNew, nEdgesNew, nVerticesNew
integer :: maxEdgesNew, maxEdges2New, TWONew, vertexDegreeNew, nVertLevelsNew
integer, allocatable, dimension(:) :: indexToCellIDNew, indexToEdgeIDNew, indexToVertexIDNew
-real, allocatable, dimension(:) :: xCellNew, yCellNew, zCellNew, latCellNew, lonCellNew
+real, allocatable, dimension(:) :: xCellNew, yCellNew, zCellNew, latCellNew, lonCellNew, meshDensityNew
real, allocatable, dimension(:) :: xEdgeNew, yEdgeNew, zEdgeNew, latEdgeNew, lonEdgeNew
real, allocatable, dimension(:) :: xVertexNew, yVertexNew, zVertexNew, latVertexNew, lonVertexNew
integer, allocatable, dimension(:) :: nEdgesOnCellNew, nEdgesOnEdgeNew, flipVerticesOnEdgeOrdering
@@ -238,6 +238,7 @@
areaCellNew, &
maxLevelCellNew, &
depthCell, &
+ temperatureNew(1,1,:), &
kiteAreasOnVertexNew )
@@ -498,7 +499,7 @@
zin = zEdgeNew(iEdge)
rlon = lonEdgeNew(iEdge)/dtr
rlat = latEdgeNew(iEdge)/dtr
- ix = nint(rlon/0.1 - 0.05) + 1
+ ix = nint(rlon/0.1 - 0.05) + nlon + 1
ix = mod(ix,nlon)+1
iy = nlat
do jcount=1,nlat
@@ -507,6 +508,9 @@
exit
endif
enddo
+ !stress = -0.1*cos(3.0*latEdgeNew(iEdge))
+ !ulon = stress
+ !ulat = 0.0
ulon = TAUX(ix,iy)
ulat = TAUY(ix,iy)
write(6,*) rlon, t_lon(ix), rlat, t_lat(iy)
@@ -543,7 +547,7 @@
if(mod(iCell,100).eq.0) write(6,*) 'l_woce t and s',iCell
rlon = lonCellNew(iCell)/dtr
rlat = latCellNew(iCell)/dtr
- ix = nint(rlon/0.1 - 0.05) + 1
+ ix = nint(rlon/0.1 - 0.05) + nlon + 1
ix = mod(ix,nlon)+1
iy = nlat
do j=1,nlat
@@ -811,6 +815,7 @@
allocate(zCell(nCells))
allocate(latCell(nCells))
allocate(lonCell(nCells))
+allocate(meshDensity(nCells))
allocate(xEdge(nEdges))
allocate(yEdge(nEdges))
allocate(zEdge(nEdges))
@@ -855,7 +860,7 @@
allocate(h(1,nVertLevels,nCells))
allocate(rho(1,nVertLevels,nCells))
-xCell=0; yCell=0; zCell=0; latCell=0; lonCell=0
+xCell=0; yCell=0; zCell=0; latCell=0; lonCell=0; meshDensity=1.0
xEdge=0; yEdge=0; zEdge=0; latEdge=0; lonEdge=0
xVertex=0; yVertex=0; zVertex=0; latVertex=0; lonVertex=0
@@ -874,6 +879,7 @@
time, &
latCell, &
lonCell, &
+ meshDensity, &
xCell, &
yCell, &
zCell, &
@@ -918,6 +924,7 @@
write(6,*) ' values from read grid, min/max'
write(6,*) ' latCell : ', minval(latCell), maxval(latCell)
write(6,*) ' lonCell : ', minval(lonCell), maxval(lonCell)
+write(6,*) ' meshDensity : ', minval(meshDensity),maxval(meshDensity)
write(6,*) ' xCell : ', minval(xCell), maxval(xCell)
write(6,*) ' yCell : ', minval(yCell), maxval(yCell)
write(6,*) ' zCell : ', minval(zCell), maxval(zCell)
@@ -996,6 +1003,7 @@
1, &
latCellNew, &
lonCellNew, &
+ meshDensityNew, &
xCellNew, &
yCellNew, &
zCellNew, &
@@ -1145,7 +1153,7 @@
do iCell=1,nCells
rlon = lonCell(iCell) / dtr
rlat = latCell(iCell) / dtr
- ix = nint(rlon*30)
+ ix = nint((rlon+180)*30) + nx
ix = mod(ix,nx)+1
iy = nint((rlat+90 )*30)
ix = max(1,ix); ix = min(nx,ix)
@@ -1190,18 +1198,6 @@
kmt(iCell) = min(kmt(iCell),kmt_neighbor_max)
enddo
-! eliminate isolated ocean cells
-do iCell=1,nCells
- flag = .true.
- if(kmt(iCell).ne.0) then
- do j=1,nEdgesOnCell(iCell)
- iCell1 = cellsOnCell(j,iCell)
- if(kmt(iCell1).ne.0) flag=.false.
- enddo
- endif
- if(flag) kmt(iCell)=0
-enddo
-
if(eliminate_inland_seas) then
call eliminateLoops(nCells,nEdges,nVertices,maxEdges,vertexDegree, &
nEdgesOnCell, cellsOnCell, verticesOnEdge, cellsOnVertex, edgesOnCell, lonCell, latCell, &
@@ -1307,6 +1303,7 @@
allocate(normalsNew(3,nEdgesNew))
allocate(latCellNew(nCellsNew))
allocate(lonCellNew(nCellsNew))
+allocate(meshDensityNew(nCellsNew))
allocate(xEdgeNew(nEdgesNew))
allocate(yEdgeNew(nEdgesNew))
allocate(zEdgeNew(nEdgesNew))
@@ -1340,7 +1337,7 @@
allocate(temperatureRestoreNew(nCellsNew))
allocate(salinityRestoreNew(nCellsNew))
-xCellNew=0; yCellNew=0; zCellNew=0; latCellNew=0; lonCellNew=0
+xCellNew=0; yCellNew=0; zCellNew=0; latCellNew=0; lonCellNew=0; meshDensityNew=1.0
xEdgeNew=0; yEdgeNew=0; zEdgeNew=0; latEdgeNew=0; lonEdgeNew=0
xVertexNew=0; yVertexNew=0; zVertexNew=0; latVertexNew=0; lonVertexNew=0
@@ -1360,6 +1357,7 @@
zCellNew(jNew)=zCell(i)
latCellNew(jNew)=latCell(i)
lonCellNew(jNew)=lonCell(i)
+ meshDensityNew(jNew)=meshDensity(i)
areaCellNew(jNew)=areaCell(i)
maxLevelCellNew(jNew) = kmt(i)
depthCell(jNew) = kmt(i)
@@ -1399,6 +1397,7 @@
deallocate(zCell)
deallocate(latCell)
deallocate(lonCell)
+deallocate(meshDensity)
deallocate(xEdge)
deallocate(yEdge)
deallocate(zEdge)
Modified: branches/ocean_projects/basin/src/module_read_netcdf.F
===================================================================
--- branches/ocean_projects/basin/src/module_read_netcdf.F 2011-10-24 20:25:34 UTC (rev 1123)
+++ branches/ocean_projects/basin/src/module_read_netcdf.F 2011-10-24 21:13:17 UTC (rev 1124)
@@ -12,6 +12,7 @@
integer :: rdDimIDvertexDegree
integer :: rdVarIDlatCell
integer :: rdVarIDlonCell
+ integer :: rdVarIDmeshDensity
integer :: rdVarIDxCell
integer :: rdVarIDyCell
integer :: rdVarIDzCell
@@ -128,6 +129,7 @@
!
nferr = nf_inq_varid(rd_ncid, 'latCell', rdVarIDlatCell)
nferr = nf_inq_varid(rd_ncid, 'lonCell', rdVarIDlonCell)
+ nferr = nf_inq_varid(rd_ncid, 'meshDensity', rdVarIDmeshDensity)
nferr = nf_inq_varid(rd_ncid, 'xCell', rdVarIDxCell)
nferr = nf_inq_varid(rd_ncid, 'yCell', rdVarIDyCell)
nferr = nf_inq_varid(rd_ncid, 'zCell', rdVarIDzCell)
@@ -175,6 +177,7 @@
time, &
latCell, &
lonCell, &
+ meshDensity, &
xCell, &
yCell, &
zCell, &
@@ -223,6 +226,7 @@
integer, intent(in) :: time
real (kind=8), dimension(:), intent(out) :: latCell
real (kind=8), dimension(:), intent(out) :: lonCell
+ real (kind=8), dimension(:), intent(out) :: meshDensity
real (kind=8), dimension(:), intent(out) :: xCell
real (kind=8), dimension(:), intent(out) :: yCell
real (kind=8), dimension(:), intent(out) :: zCell
@@ -262,12 +266,16 @@
real (kind=8), dimension(:,:,:), intent(out) :: u
real (kind=8), dimension(:,:,:), intent(out) :: v
real (kind=8), dimension(:,:,:), intent(out) :: h
+
+ logical :: meshDensityPresent
integer :: nferr
integer, dimension(1) :: start1, count1
integer, dimension(2) :: start2, count2
integer, dimension(3) :: start3, count3
integer, dimension(4) :: start4, count4
+
+ meshDensityPresent = .false.
start1(1) = 1
@@ -287,6 +295,17 @@
count1( 1) = rdLocalnCells
count1( 1) = rdLocalnCells
nferr = nf_get_vara_double(rd_ncid, rdVarIDlonCell, start1, count1, lonCell)
+
+ start1(1) = 1
+ count1( 1) = rdLocalnCells
+ count1( 1) = rdLocalnCells
+ nferr = nf_inq_varid(rd_ncid, 'meshDensity', rdVarIDmeshDensity)
+ if(nferr.eq.0) then
+ nferr = nf_get_vara_double(rd_ncid, rdVarIDmeshDensity, start1, count1, meshDensity)
+ else
+ meshDensity=1.0
+ write(6,*) ' mesh density not present ', nferr, rdVarIDmeshDensity
+ endif
start1(1) = 1
count1( 1) = rdLocalnCells
Modified: branches/ocean_projects/basin/src/module_write_netcdf.F
===================================================================
--- branches/ocean_projects/basin/src/module_write_netcdf.F 2011-10-24 20:25:34 UTC (rev 1123)
+++ branches/ocean_projects/basin/src/module_write_netcdf.F 2011-10-24 21:13:17 UTC (rev 1124)
@@ -12,6 +12,7 @@
integer :: wrDimIDnVertLevels
integer :: wrVarIDlatCell
integer :: wrVarIDlonCell
+ integer :: wrVarIDmeshDensity
integer :: wrVarIDxCell
integer :: wrVarIDyCell
integer :: wrVarIDzCell
@@ -130,6 +131,8 @@
dimlist( 1) = wrDimIDnCells
nferr = nf_def_var(wr_ncid, 'lonCell', NF_DOUBLE, 1, dimlist, wrVarIDlonCell)
dimlist( 1) = wrDimIDnCells
+ nferr = nf_def_var(wr_ncid, 'meshDensity', NF_DOUBLE, 1, dimlist, wrVarIDmeshDensity)
+ dimlist( 1) = wrDimIDnCells
nferr = nf_def_var(wr_ncid, 'xCell', NF_DOUBLE, 1, dimlist, wrVarIDxCell)
dimlist( 1) = wrDimIDnCells
nferr = nf_def_var(wr_ncid, 'yCell', NF_DOUBLE, 1, dimlist, wrVarIDyCell)
@@ -269,6 +272,7 @@
time, &
latCell, &
lonCell, &
+ meshDensity, &
xCell, &
yCell, &
zCell, &
@@ -327,6 +331,7 @@
integer, intent(in) :: time
real (kind=8), dimension(:), intent(in) :: latCell
real (kind=8), dimension(:), intent(in) :: lonCell
+ real (kind=8), dimension(:), intent(in) :: meshDensity
real (kind=8), dimension(:), intent(in) :: xCell
real (kind=8), dimension(:), intent(in) :: yCell
real (kind=8), dimension(:), intent(in) :: zCell
@@ -405,6 +410,10 @@
start1(1) = 1
count1( 1) = wrLocalnCells
nferr = nf_put_vara_double(wr_ncid, wrVarIDlonCell, start1, count1, lonCell)
+
+ start1(1) = 1
+ count1( 1) = wrLocalnCells
+ nferr = nf_put_vara_double(wr_ncid, wrVarIDmeshDensity, start1, count1, meshDensity)
start1(1) = 1
count1( 1) = wrLocalnCells
Modified: branches/ocean_projects/basin/src/utilities.F
===================================================================
--- branches/ocean_projects/basin/src/utilities.F 2011-10-24 20:25:34 UTC (rev 1123)
+++ branches/ocean_projects/basin/src/utilities.F 2011-10-24 21:13:17 UTC (rev 1124)
@@ -25,6 +25,7 @@
areaCell, &
maxLevelCell, &
depthCell, &
+ SST, &
kiteAreasOnVertex )
implicit none
@@ -47,13 +48,13 @@
integer, dimension(vertexDegree, nVertices), intent(in) :: cellsOnVertex
integer, dimension(nCells), intent(in) :: maxLevelCell
real (kind=8), dimension(nCells), intent(in) :: areaCell
- real (kind=8), dimension(nCells), intent(in) :: depthCell
+ real (kind=8), dimension(nCells), intent(in) :: depthCell, SST
real (kind=8), dimension(vertexDegree,nVertices), intent(in) :: kiteAreasOnVertex
character(len=80) :: a, b, c, d, e, f
integer :: i, j, k, nVerticesTotal, iEdge, iLoop, iFace, Vert(4), Edge(4), iVertex, i1, i2, jp1
integer :: nKitesTotal, iCell, iEdge1, iEdge2, iVertex11, iVertex12, iVertex21, iVertex22, ksave
- real (kind=8) :: x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, xscale, work(nCells)
+ real (kind=8) :: x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, xscale, work(nCells), work1(nCells)
real (kind=8) :: xv, yv, zv, xc, yc, zc, dist
logical (kind=8) :: eflag
@@ -173,8 +174,10 @@
close(1)
- work = depthCell
+ work1 = depthCell
+ work = SST
+ open(unit= 9,file='dx/ocean.depth.data',form='formatted',status='unknown')
open(unit=10,file='dx/ocean.area.data',form='formatted',status='unknown')
open(unit=11,file='dx/ocean.face.data',form='formatted',status='unknown')
open(unit=12,file='dx/ocean.loop.data',form='formatted',status='unknown')
@@ -184,6 +187,7 @@
iLoop = 0
iEdge = 0
do i=1,nCells
+ write(9,20) work1(i)
write(10,20) work(i)
write(11,21) i-1
write(12,21) iLoop
@@ -234,6 +238,7 @@
22 format(3e20.10)
23 format(3i8, 3e20.10)
+ close(9)
close(10)
close(11)
close(12)
</font>
</pre>