<p><b>duda</b> 2012-02-06 10:35:05 -0700 (Mon, 06 Feb 2012)</p><p>BRANCH COMMIT<br>
<br>
Several small fixes for the real-data initialization case:<br>
<br>
- use sst_prefix rather than met_prefix for the source of SST and SEAICE fields<br>
to go in the surface update files<br>
<br>
- correct grid spacing in polar stereographic projection (no need to convert km to m)<br>
<br>
- fix logic so that, when intermediate files have both SST and SEAICE, both fields<br>
will be interpolated<br>
<br>
<br>
M src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F<br>
</p><hr noshade><pre><font color="gray">Modified: branches/atmos_physics/src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F
===================================================================
--- branches/atmos_physics/src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F        2012-02-03 21:33:51 UTC (rev 1465)
+++ branches/atmos_physics/src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F        2012-02-06 17:35:05 UTC (rev 1466)
@@ -3883,7 +3883,7 @@
if (field % iproj == PROJ_PS) then
call map_set(PROJ_PS, proj, &
- dx = real(field % dx * 1000.0,RKIND), &
+ dx = real(field % dx,RKIND), &
truelat1 = real(field % truelat1,RKIND), &
stdlon = real(field % xlonc,RKIND), &
knowni = real(field % nx / 2.0,RKIND), &
@@ -4346,12 +4346,12 @@
curr_time = mpas_get_clock_time(fg_clock, MPAS_NOW)
do while (curr_time <= stop_time)
call mpas_get_time(curr_time, dateTimeString=timeString)
- write(0,*) 'Processing ',trim(config_met_prefix)//':'//timeString(1:13)
+ write(0,*) 'Processing ',trim(config_sst_prefix)//':'//timeString(1:13)
! Open intermediate file
- call read_met_init(trim(config_met_prefix), .false., timeString(1:13), istatus)
+ call read_met_init(trim(config_sst_prefix), .false., timeString(1:13), istatus)
if (istatus /= 0) then
- write(0,*) 'Error reading ',trim(config_met_prefix)//':'//timeString(1:13)
+ write(0,*) 'Error reading ',trim(config_sst_prefix)//':'//timeString(1:13)
exit
end if
@@ -4389,6 +4389,15 @@
lat1 = real(field % startlat,RKIND), &
lon1 = real(field % startlon,RKIND))
! nxmax = nint(360.0 / field % deltalon), &
+ else if (field % iproj == PROJ_PS) then
+ call map_set(PROJ_PS, proj, &
+ dx = real(field % dx,RKIND), &
+ truelat1 = real(field % truelat1,RKIND), &
+ stdlon = real(field % xlonc,RKIND), &
+ knowni = real(field % nx / 2.0,RKIND), &
+ knownj = real(field % ny / 2.0,RKIND), &
+ lat1 = real(field % startlat,RKIND), &
+ lon1 = real(field % startlon,RKIND))
end if
! Interpolate SST/SKINTEMP field to each MPAS grid cell
@@ -4413,10 +4422,8 @@
deallocate(slab_r8)
deallocate(field % slab)
- exit
- end if
- if (index(field % field, 'SEAICE') /= 0) then
+ else if (index(field % field, 'SEAICE') /= 0) then
! Interpolation routines use real(kind=RKIND), so copy from default real array
allocate(slab_r8(field % nx, field % ny))
@@ -4446,6 +4453,15 @@
lat1 = real(field % startlat,RKIND), &
lon1 = real(field % startlon,RKIND))
! nxmax = nint(360.0 / field % deltalon), &
+ else if (field % iproj == PROJ_PS) then
+ call map_set(PROJ_PS, proj, &
+ dx = real(field % dx,RKIND), &
+ truelat1 = real(field % truelat1,RKIND), &
+ stdlon = real(field % xlonc,RKIND), &
+ knowni = real(field % nx / 2.0,RKIND), &
+ knownj = real(field % ny / 2.0,RKIND), &
+ lat1 = real(field % startlat,RKIND), &
+ lon1 = real(field % startlon,RKIND))
end if
! Interpolate SEAICE/SKINTEMP field to each MPAS grid cell
@@ -4466,15 +4482,18 @@
call latlon_to_ij(proj, lat, lon, x, y)
end if
fg % xice % array(iCell) = interp_sequence(x, y, 1, slab_r8, 1, field % nx, 1, field % ny, 1, 1, -1.e30_RKIND, interp_list, 1)
+ if (fg % xice % array(iCell) == -1.e30_RKIND) fg % xice % array(iCell) = 0.0_RKIND
end do
deallocate(slab_r8)
deallocate(field % slab)
- exit
+
+ else
+
+ deallocate(field % slab)
end if
- deallocate(field % slab)
call read_next_met_field(field, istatus)
end do
</font>
</pre>