1c1
< !  Create an initial data set for the WRF model based on real data.  This
---
> !i  Create an initial data set for the WRF model based on real data.  This
5d4
<    USE module_big_step_utilities_em
11d9
<    USE module_bc
14a13
>    USE esmf_mod
21c20
<     REAL    :: time , bdyfrq
---
>    REAL    :: time , bdyfrq
23c22
<    INTEGER :: loop , levels_to_process
---
>    INTEGER :: loop , levels_to_process , debug_level
26c25,26
<    TYPE(domain) , POINTER :: keep_grid, grid_ptr, null_domain, grid
---
>    TYPE(domain) , POINTER :: null_domain
>    TYPE(domain) , POINTER :: grid
31d30
<    INTEGER :: id1 , id , ierr
42c41,43
<    INTEGER :: i , j , k , idts
---
>    INTEGER :: ips , ipe , jps , jpe , kps , kpe
>    INTEGER :: ijds , ijde , spec_bdy_width
>    INTEGER :: i , j , k , idts, rc
46,54d46
<    REAL , DIMENSION(:,:,:) , ALLOCATABLE :: ubdy3dtemp1 , vbdy3dtemp1 , tbdy3dtemp1 , pbdy3dtemp1 , qbdy3dtemp1
<    REAL , DIMENSION(:,:  ) , ALLOCATABLE :: mbdy2dtemp1
<    REAL , DIMENSION(:,:,:) , ALLOCATABLE :: ubdy3dtemp2 , vbdy3dtemp2 , tbdy3dtemp2 , pbdy3dtemp2 , qbdy3dtemp2
<    REAL , DIMENSION(:,:  ) , ALLOCATABLE :: mbdy2dtemp2
< 
<    CHARACTER(LEN=24) :: previous_date , this_date , next_date
<    CHARACTER(LEN=19) :: start_date_char , end_date_char , current_date_char , next_date_char
<    CHARACTER(LEN= 4) :: loop_char
< 
58a51,58
>    TYPE(ESMF_TimeInterval) :: time_interval
> real::t1,t2
>    INTERFACE
>      SUBROUTINE Set_Timekeeping( grid )
>       USE module_domain
>       TYPE(domain), POINTER :: grid
>      END SUBROUTINE Set_Timekeeping
>    END INTERFACE
60,64d59
<    CHARACTER (LEN=80) :: inpname , bdyname
< 
< 
<    !  Get the NAMELIST data for input.
< 
67c62
<    program_name = "REAL_EM V1.3 PREPROCESSOR"
---
>    program_name = "REAL_EM V1.4 PREPROCESSOR"
69d63
< 
74c68,71
<    CALL init_modules_em
---
>    !  Initialize the modules used by the WRF system.  Many of the CALLs made from the
>    !  init_modules routine are NO-OPs.  Typical initializations are: the size of a
>    !  REAL, setting the file handles to a pre-use value, defining moisture and
>    !  chemistry indices, etc.
75a73,77
>    CALL       wrf_debug ( 100 , 'real_em: calling init_modules ' )
>    CALL init_modules
> 
>    !  The configuration switches mostly come from the NAMELIST input.
> 
78,81c80
<      CALL initial_config
<      CALL get_config_as_buffer( configbuf, configbuflen, nbytes )
<      CALL wrf_dm_bcast_bytes( configbuf, nbytes )
<      CALL set_config_as_buffer( configbuf, configbuflen )
---
>       CALL initial_config
82a82,84
>    CALL get_config_as_buffer( configbuf, configbuflen, nbytes )
>    CALL wrf_dm_bcast_bytes( configbuf, nbytes )
>    CALL set_config_as_buffer( configbuf, configbuflen )
87a90,92
>    CALL get_debug_level ( debug_level )
>    CALL set_wrf_debug_level ( debug_level )
> 
90c95
<    !  An available simple timer from the timing module.
---
>    !  Allocate the space for the mother of all domains.
92a98
>    CALL       wrf_debug ( 100 , 'real_em: calling alloc_and_configure_domain ' )
94d99
<                                      local_time = 0           , &
98d102
< print *,'ok'
100d103
<    CALL set_scalar_indices_from_config ( head_grid%id , idum1, idum2 )
102a106,113
>    CALL Set_Timekeeping ( grid )
>    CALL ESMF_TimeIntervalSet ( time_interval , S=model_config_rec%interval_seconds, rc=rc )
>    CALL ESMF_ClockSetTimeStep ( grid%domain_clock , time_interval , rc=rc )
> 
>    CALL       wrf_debug ( 100 , 'real_em: calling set_scalar_indices_from_config ' )
>    CALL set_scalar_indices_from_config ( grid%id , idum1, idum2 )
> 
>    CALL       wrf_debug ( 100 , 'real_em: calling model_to_grid_config_rec ' )
105,110c116
< print *,'start date=',model_config_rec%start_year(grid%id),model_config_rec%start_month(grid%id),&
< model_config_rec%start_day(grid%id),model_config_rec%start_hour(grid%id)
< print *,'end   date=',model_config_rec%end_year(grid%id),model_config_rec%end_month(grid%id),&
< model_config_rec%end_day(grid%id),model_config_rec%end_hour(grid%id)
< print *,'interval  =',model_config_rec%interval_seconds
< print *,'init_typ  =',model_config_rec%real_data_init_type
---
>    !  Initialize the WRF IO: open files, init file handles, etc.
112c118,119
<    !  Figure out the starting and ending dates in a character format.
---
>    CALL       wrf_debug ( 100 , 'real_em: calling init_wrfio' )
>    CALL init_wrfio
114,119c121,122
<    start_year   = model_config_rec%start_year  (grid%id)
<    start_month  = model_config_rec%start_month (grid%id)
<    start_day    = model_config_rec%start_day   (grid%id)
<    start_hour   = model_config_rec%start_hour  (grid%id)
<    start_minute = model_config_rec%start_minute(grid%id)
<    start_second = model_config_rec%start_second(grid%id)
---
>    !  Some of the configuration values may have been modified from the initial READ
>    !  of the NAMELIST, so we re-broadcast the configuration records.
121,126c124,129
<      end_year   = model_config_rec%  end_year  (grid%id)
<      end_month  = model_config_rec%  end_month (grid%id)
<      end_day    = model_config_rec%  end_day   (grid%id)
<      end_hour   = model_config_rec%  end_hour  (grid%id)
<      end_minute = model_config_rec%  end_minute(grid%id)
<      end_second = model_config_rec%  end_second(grid%id)
---
> #ifdef DM_PARALLEL
>    CALL       wrf_debug ( 100 , 'real_em: re-broadcast the configuration records' )
>    CALL get_config_as_buffer( configbuf, configbuflen, nbytes )
>    CALL wrf_dm_bcast_bytes( configbuf, nbytes )
>    CALL set_config_as_buffer( configbuf, configbuflen )
> #endif
128,129c131
<    interval_seconds    = model_config_rec%interval_seconds
<    real_data_init_type = model_config_rec%real_data_init_type
---
>    !   No looping in this layer.  
130a133,319
>    CALL       wrf_debug ( 100 , 'calling med_sidata_input' )
>    CALL med_sidata_input ( grid , config_flags )
>    CALL       wrf_debug ( 100 , 'backfrom med_sidata_input' )
> 
>    !  We are done.
> 
>    CALL       wrf_debug (   0 , 'real_em: SUCCESS COMPLETE REAL_EM INIT' )
> 
> #ifdef DM_PARALLEL
>     CALL wrf_dm_shutdown
> #endif
> 
> END PROGRAM real_data
> 
> SUBROUTINE med_sidata_input ( grid , config_flags )
>   ! Driver layer
>    USE module_domain
>    USE module_io_domain
>   ! Model layer
>    USE module_configure
>    USE module_bc_time_utilities
>    USE module_initialize
>    USE module_optional_si_input
> 
>    USE module_date_time
> 
>    IMPLICIT NONE
> 
> 
>   ! Interface 
>    INTERFACE
>      SUBROUTINE start_domain ( grid )  ! comes from module_start in appropriate dyn_ directory
>        USE module_domain
>        TYPE (domain) grid
>      END SUBROUTINE start_domain
>    END INTERFACE
> 
>   ! Arguments
>    TYPE(domain)                :: grid
>    TYPE (grid_config_rec_type) :: config_flags
>   ! Local
>    INTEGER                :: time_step_begin_restart
>    INTEGER                :: idsi , ierr , myproc
>    CHARACTER (LEN=80)      :: si_inpname
>    CHARACTER (LEN=80)      :: message
> 
>    CHARACTER(LEN=19) :: start_date_char , end_date_char , current_date_char , next_date_char
> 
>    INTEGER :: time_loop_max , loop, rc
>    INTEGER :: julyr , julday 
>    REAL :: gmt
> real::t1,t2
> 
>    grid%input_from_file = .true.
>    grid%input_from_file = .false.
> 
>    CALL compute_si_start_and_end ( model_config_rec%start_year  (grid%id) , &
>                                    model_config_rec%start_month (grid%id) , &
>                                    model_config_rec%start_day   (grid%id) , &
>                                    model_config_rec%start_hour  (grid%id) , &
>                                    model_config_rec%start_minute(grid%id) , &
>                                    model_config_rec%start_second(grid%id) , &
>                                    model_config_rec%  end_year  (grid%id) , & 
>                                    model_config_rec%  end_month (grid%id) , &
>                                    model_config_rec%  end_day   (grid%id) , &
>                                    model_config_rec%  end_hour  (grid%id) , &
>                                    model_config_rec%  end_minute(grid%id) , &
>                                    model_config_rec%  end_second(grid%id) , &
>                                    model_config_rec%interval_seconds      , &
>                                    model_config_rec%real_data_init_type   , &
>                                    start_date_char , end_date_char , time_loop_max )
> 
>    !  Here we define the initial time to process, for later use by the code.
>    
>    print*,'line 207: start_date_char=',start_date_char,'  current_date_char=',current_date_char
>    current_date_char = start_date_char
>    start_date = start_date_char // '.0000'
>    current_date = start_date
>    print*,'line 211: start_date=',start_date,'  current_date=',current_date
> 
>    CALL set_bdyfrq ( grid%id , REAL(model_config_rec%interval_seconds) )
> 
>    !!!!!!!  Loop over each time period to process.
> 
>    DO loop = 1 , time_loop_max
> 
>       print *,'-----------------------------------------------------------------------------'
>       print *,' '
>       print '(A,A,A,I2,A,I2)' , ' Current date being processed: ',current_date, ', which is loop #',loop,' out of ',time_loop_max
> 
>       !  After current_date has been set, fill in the julgmt stuff.
> 
>       CALL geth_julgmt ( config_flags%julyr , config_flags%julday , config_flags%gmt )
> 
> 	print *,'configflags%julyr, %julday, %gmt:',config_flags%julyr, config_flags%julday, config_flags%gmt
>       !  Now that the specific Julian info is available, save these in the model config record.
> 
>       CALL set_gmt (grid%id, config_flags%gmt)
>       CALL set_julyr (grid%id, config_flags%julyr)
>       CALL set_julday (grid%id, config_flags%julday)
> 
>       !  Open the wrfinput file.
> 
>       IF ( grid%dyn_opt .EQ. dyn_em ) THEN
>          CALL       wrf_debug ( 100 , 'med_sidata_input: calling open_r_dataset for wrf_real_input_em' )
>          CALL construct_filename2( si_inpname , 'wrf_real_input_em' , grid%id , 2 , current_date_char )
>       END IF
> 	write(0,*) 'call open_r_dataset'
>       CALL open_r_dataset ( idsi, TRIM(si_inpname) , grid , config_flags , "DATASET=INPUT", ierr )
> 	write(0,*) 'return open_r_dataset'
>       IF ( ( ierr .NE. 0 ) .AND. ( grid%dyn_opt .EQ. dyn_em ) ) THEN
>         CALL wrf_error_fatal( 'real: error opening wrf_real_input_em for reading' )
>       ENDIF
> 
> 	write(0,*) 'here 1'
> 
>       !  Input data.
> 
>       CALL       wrf_debug ( 100 , 'med_sidata_input: calling input_aux_model_input1_wrf' )
>       CALL input_aux_model_input1_wrf ( idsi ,   grid , config_flags , ierr )
> 	write(0,*) 'here 2'
> 
>       !  Possible optional SI input.  This sets flags used by init_domain.
>       IF ( loop .EQ. 1 ) THEN
>          CALL       wrf_debug ( 100 , 'med_sidata_input: calling init_module_optional_si_input' )
> 	write(0,*) 'here 3'
>          CALL init_module_optional_si_input ( grid , config_flags )
> 	write(0,*) 'here 4'
>       END IF
>       CALL       wrf_debug ( 100 , 'med_sidata_input: calling optional_si_input' )
> 	write(0,*) 'here 5'
>       CALL  optional_si_input ( grid , idsi )
> 	write(0,*) 'here 6'
> 
>       !  Initialize the mother domain for this time period with input data.
> 
>       CALL       wrf_debug ( 100 , 'med_sidata_input: calling init_domain' )
>       grid%input_from_file = .true.
> 	write(0,*) 'here 7'
>       CALL init_domain (  grid )
> 	write(0,*) 'here 8'
>       CALL model_to_grid_config_rec ( grid%id , model_config_rec , config_flags )
> 	write(0,*) 'here 9'
> 
>       !  Close this file that is output from the SI and input to this pre-proc.
> 
>       CALL       wrf_debug ( 100 , 'med_sidata_input: back from init_domain' )
> 	write(0,*) 'here 10'
>       CALL close_dataset ( idsi , config_flags , "DATASET=INPUT" )
> 
> !     CALL start_domain ( grid )
> 
> 	write(0,*) 'call assemble_output'
> 
>       CALL assemble_output ( grid , config_flags , loop , time_loop_max )
> 	write(0,*) 'return assemble_output'
> 
>       !  Here we define the next time that we are going to process.
> 
>       CALL geth_newdate ( current_date_char , start_date_char , loop * model_config_rec%interval_seconds )
>       current_date =  current_date_char // '.0000'
>       CALL atotime( current_date(1:19), grid%current_time )
>       CALL ESMF_ClockSetCurrTime(grid%domain_clock, grid%current_time, rc)
> if (loop.ne.time_loop_max )then
> print *,'set esmf clock to ',current_date(1:19)
> endif
> 
>    END DO
> END SUBROUTINE med_sidata_input
> 
> SUBROUTINE compute_si_start_and_end (  &
>    start_year , start_month , start_day , start_hour , start_minute , start_second , &
>      end_year ,   end_month ,   end_day ,   end_hour ,   end_minute ,   end_second , &
>    interval_seconds , real_data_init_type , &
>    start_date_char , end_date_char , time_loop_max )
> 
>    USE module_date_time
> 
>    IMPLICIT NONE
> 
>    INTEGER :: start_year , start_month , start_day , start_hour , start_minute , start_second
>    INTEGER ::   end_year ,   end_month ,   end_day ,   end_hour ,   end_minute ,   end_second
>    INTEGER :: interval_seconds , real_data_init_type
>    INTEGER :: time_loop_max , time_loop
> 
>    CHARACTER(LEN=19) :: current_date_char , start_date_char , end_date_char , next_date_char
> 
135a325,326
> !  start_date = start_date_char // '.0000'
> 
144c335
<          time_loop = time_loop + 1 
---
>          time_loop = time_loop + 1
148c339
<          time_loop = time_loop + 1 
---
>          time_loop = time_loop + 1
158a350
> END SUBROUTINE compute_si_start_and_end
160c352
<    !  Here we define the initial time to process, for later use by the code.
---
> SUBROUTINE assemble_output ( grid , config_flags , loop , time_loop_max )
162,168c354,360
<    current_date_char = start_date_char
<    start_date = start_date_char // '.0000'
<    current_date = start_date
<    bdyfrq = interval_seconds
<    CALL set_bdyfrq ( grid%id , bdyfrq )
<   
<    !  Loop over each time period to process.
---
>    USE module_big_step_utilities_em
>    USE module_domain
>    USE module_io_domain
>    USE module_configure
>    USE module_date_time
>    USE module_bc
>    IMPLICIT NONE
170c362,364
<    DO loop = 1 , time_loop_max
---
>    TYPE(domain)                 :: grid
>    TYPE (grid_config_rec_type)  :: config_flags
>    INTEGER , INTENT(IN)         :: loop , time_loop_max
172c366,370
<       !  Initialize the mother domain for this time period.
---
>    INTEGER :: ids , ide , jds , jde , kds , kde
>    INTEGER :: ims , ime , jms , jme , kms , kme
>    INTEGER :: ips , ipe , jps , jpe , kps , kpe
>    INTEGER :: ijds , ijde , spec_bdy_width
>    INTEGER :: i , j , k , idts
174,180c372,377
<       grid%input_from_file = .false.
<       CALL init_domain (  grid )
<       CALL model_to_grid_config_rec ( grid%id , model_config_rec , config_flags )
<    
<       IF ( loop .EQ. 1 ) THEN
<          CALL init_wrfio
<       END IF
---
>    INTEGER :: id1 , interval_seconds , ierr, rc
>    INTEGER , SAVE :: id 
>    CHARACTER (LEN=80) :: inpname , bdyname
>    CHARACTER(LEN= 4) :: loop_char
> character *19 :: temp19
> character *24 :: temp24 , temp24b
182,187c379,383
<       !  There are different things that have to happen during different time loops.
<       !  During EVERY time loop, EXCEPT the last, we output the starting values for the
<       !  boundary data.  During EVERY time loop, EXCEPT the first, we output the 
<       !  boundary tendency.
<    
<       IF ( loop .EQ. 1 ) THEN
---
>    REAL , DIMENSION(:,:,:) , ALLOCATABLE , SAVE :: ubdy3dtemp1 , vbdy3dtemp1 , tbdy3dtemp1 , pbdy3dtemp1 , qbdy3dtemp1
>    REAL , DIMENSION(:,:  ) , ALLOCATABLE , SAVE :: mbdy2dtemp1
>    REAL , DIMENSION(:,:,:) , ALLOCATABLE , SAVE :: ubdy3dtemp2 , vbdy3dtemp2 , tbdy3dtemp2 , pbdy3dtemp2 , qbdy3dtemp2
>    REAL , DIMENSION(:,:  ) , ALLOCATABLE , SAVE :: mbdy2dtemp2
> real::t1,t2
189c385
<          !  Open the wrfinput file.
---
>    !  Various sizes that we need to be concerned about.
191,214c387,392
<          CALL construct_filename1( inpname , 'wrfinput' , grid%id , 2 )
<          CALL open_w_dataset ( id1, TRIM(inpname) , head_grid , config_flags , output_model_input , "DATASET=INPUT", ierr )
<          IF ( ierr .NE. 0 ) THEN
<            CALL wrf_error_fatal( 'real: error opening wrfinput for writing' )
<          ENDIF
<    
<          !  Various sizes that we need to be concerned about.
<    
<          ids = grid%sd31 
<          ide = grid%ed31
<          kds = grid%sd32 
<          kde = grid%ed32
<          jds = grid%sd33 
<          jde = grid%ed33
<       
<          ims = grid%sm31 
<          ime = grid%em31
<          kms = grid%sm32 
<          kme = grid%em32
<          jms = grid%sm33 
<          jme = grid%em33
<    
<          print *, ids , ide , jds , jde , kds , kde
<          print *, ims , ime , jms , jme , kms , kme
---
>    ids = grid%sd31
>    ide = grid%ed31
>    kds = grid%sd32
>    kde = grid%ed32
>    jds = grid%sd33
>    jde = grid%ed33
216,217c394,399
<          !  This is the space needed to save the current 3d data for use in computing
<          !  the lateral boundary tendencies.
---
>    ims = grid%sm31
>    ime = grid%em31
>    kms = grid%sm32
>    kme = grid%em32
>    jms = grid%sm33
>    jme = grid%em33
219,230c401,406
<          ALLOCATE ( ubdy3dtemp1(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( vbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( tbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( pbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( qbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( mbdy2dtemp1(ims:ime,        jms:jme) )
<          ALLOCATE ( ubdy3dtemp2(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( vbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( tbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( pbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( qbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
<          ALLOCATE ( mbdy2dtemp2(ims:ime,        jms:jme) )
---
>    ips = grid%sp31
>    ipe = grid%ep31
>    kps = grid%sp32
>    kpe = grid%ep32
>    jps = grid%sp33
>    jpe = grid%ep33
232c408,409
<          !  Actually do work here - output the data!
---
>    ijds = MIN ( ids , jds )
>    ijde = MAX ( ide , jde )
234,237c411
<          CALL calc_current_date ( grid%id , 0. )
<          head_grid%write_metadata = .true.
<          CALL output_model_input ( id1, head_grid , config_flags , ierr )
<          CALL close_dataset ( id1 , config_flags , "DATASET=INPUT" )
---
>    !  Boundary width, scalar value.
239,252c413,414
<          !  We need to save the 3d data to compute a difference during the next loop.  Couple the
<          !  3d fields with total mu (mub + mu_1) and the stagger-specific map scale factor.
<       
<          CALL couple ( grid%em_mu_1 , grid%em_mub , ubdy3dtemp1 , grid%em_u_1                 , 'u' , grid%msfu , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , vbdy3dtemp1 , grid%em_v_1                 , 'v' , grid%msfv , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , tbdy3dtemp1 , grid%em_t_1                 , 't' , grid%msft , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , pbdy3dtemp1 , grid%em_ph_1                , 'w' , grid%msft , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , qbdy3dtemp1 , grid%moist_1(:,:,:,P_QV) , 't' , grid%msft , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          mbdy2dtemp1 = grid%em_mu_1
---
>    spec_bdy_width = model_config_rec%spec_bdy_width
>    interval_seconds = model_config_rec%interval_seconds
254,285c416
<          !  There are 2 components to the lateral boundaries.  First, there is the starting
<          !  point of this time period - just the outer few rows and columns.
<    
<          CALL stuff_bdy     ( ubdy3dtemp1 , grid%em_u_b     , 'U' , ids , ide , jds , jde , kds , kde , &
<                                                                 ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdy     ( vbdy3dtemp1 , grid%em_v_b     , 'V' , ids , ide , jds , jde , kds , kde , &
<                                                                 ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdy     ( tbdy3dtemp1 , grid%em_t_b     , 'T' , ids , ide , jds , jde , kds , kde , &
<                                                                 ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdy     ( pbdy3dtemp1 , grid%em_ph_b    , 'W' , ids , ide , jds , jde , kds , kde , &
<                                                                 ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdy     ( qbdy3dtemp1 , grid%em_rqv_b   , 'T' , ids , ide , jds , jde , kds , kde , &
<                                                                 ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdy2    ( mbdy2dtemp1 , grid%em_mu_b    , 'T' , ids , ide , jds , jde , kds , kde , &
<                                                                 ims , ime , jms , jme , kms , kme ) 
<       ELSE IF ( loop .GT. 1 ) THEN
<       
<          CALL couple ( grid%em_mu_1 , grid%em_mub , ubdy3dtemp2 , grid%em_u_1                 , 'u' , grid%msfu , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , vbdy3dtemp2 , grid%em_v_1                 , 'v' , grid%msfv , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , tbdy3dtemp2 , grid%em_t_1                 , 't' , grid%msft , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , pbdy3dtemp2 , grid%em_ph_1                , 'w' , grid%msft , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          CALL couple ( grid%em_mu_1 , grid%em_mub , qbdy3dtemp2 , grid%moist_1(:,:,:,P_QV) , 't' , grid%msft , &
<                        ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ids, ide, jds, jde, kds, kde )
<          mbdy2dtemp2 = grid%em_mu_1
<    
<          !  During all of the loops after the first loop, we first compute the boundary
<          !  tendencies with the current data values and the previously save information
<          !  stored in the *bdy3dtemp1 arrays.
---
> print *,'spec_bdy_width, interval_seconds, loop=',spec_bdy_width , interval_seconds,loop
287,306c418
<          CALL stuff_bdytend ( ubdy3dtemp2 , ubdy3dtemp1 , REAL(interval_seconds) , grid%em_u_bt  , 'U' , &
<                                                                ids , ide , jds , jde , kds , kde , &
<                                                                ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdytend ( vbdy3dtemp2 , vbdy3dtemp1 , REAL(interval_seconds) , grid%em_v_bt  , 'V' , &
<                                                                ids , ide , jds , jde , kds , kde , &
<                                                                ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdytend ( tbdy3dtemp2 , tbdy3dtemp1 , REAL(interval_seconds) , grid%em_t_bt  , 'T' , &
<                                                                ids , ide , jds , jde , kds , kde , &
<                                                                ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdytend ( pbdy3dtemp2 , pbdy3dtemp1 , REAL(interval_seconds) , grid%em_ph_bt  , 'W' , &
<                                                                ids , ide , jds , jde , kds , kde , &
<                                                                ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdytend ( qbdy3dtemp2 , qbdy3dtemp1 , REAL(interval_seconds) , grid%em_rqv_bt , 'T' , &
<                                                                ids , ide , jds , jde , kds , kde , &
<                                                                ims , ime , jms , jme , kms , kme ) 
<          CALL stuff_bdytend2( mbdy2dtemp2 , mbdy2dtemp1 , REAL(interval_seconds) , grid%em_mu_bt  , 'T' , &
<                                                                ids , ide , jds , jde , kds , kde , &
<                                                                ims , ime , jms , jme , kms , kme ) 
<                                                                
<          !  Open the boundary file.
---
>    IF ( loop .EQ. 1 ) THEN
308,317c420,421
<          IF ( loop .eq. 2 ) THEN
<             CALL construct_filename1( bdyname , 'wrfbdy' , grid%id , 2 )
<             CALL open_w_dataset ( id, TRIM(bdyname) , head_grid , config_flags , output_boundary , "DATASET=BOUNDARY", ierr )
<             IF ( ierr .NE. 0 ) THEN
<                   CALL wrf_error_fatal( 'real: error opening wrfbdy for writing' )
<             ENDIF
<             head_grid%write_metadata = .true.
<          ELSE
<             head_grid%write_metadata = .false.
<          END IF
---
>       !  This is the space needed to save the current 3d data for use in computing
>       !  the lateral boundary tendencies.
319c423,434
<          !  Both pieces of the boundary data are now available to be written.
---
>       ALLOCATE ( ubdy3dtemp1(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( vbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( tbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( pbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( qbdy3dtemp1(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( mbdy2dtemp1(ims:ime,        jms:jme) )
>       ALLOCATE ( ubdy3dtemp2(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( vbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( tbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( pbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( qbdy3dtemp2(ims:ime,kms:kme,jms:jme) )
>       ALLOCATE ( mbdy2dtemp2(ims:ime,        jms:jme) )
321c436
<          CALL output_boundary ( id, head_grid , config_flags , ierr )
---
>       !  Open the wrfinput file.  From this program, this is an *output* file.
323,324c438,446
<          !  OK, for all of the loops, we output the initialzation data, which would allow us to 
<          !  start the model at any of the available analysis time periods.
---
>       CALL construct_filename1( inpname , 'wrfinput' , grid%id , 2 )
>       CALL open_w_dataset ( id1, TRIM(inpname) , grid , config_flags , output_model_input , "DATASET=INPUT", ierr )
>       IF ( ierr .NE. 0 ) THEN
>          CALL wrf_error_fatal( 'real: error opening wrfinput for writing' )
>       ENDIF
> !     CALL calc_current_date ( grid%id , 0. )
>       grid%write_metadata = .true.
>       CALL output_model_input ( id1, grid , config_flags , ierr )
>       CALL close_dataset ( id1 , config_flags , "DATASET=INPUT" )
326,331c448,449
<          WRITE ( loop_char , FMT = '(I4.4)' ) loop
< !        CALL open_w_dataset ( id1, 'wrfinput'//loop_char , head_grid , config_flags , output_model_input , "DATASET=INPUT", ierr )
< !        IF ( ierr .NE. 0 ) THEN
< !          CALL wrf_error_fatal( 'real: error opening wrfinput'//loop_char//' for writing' )
< !        ENDIF
<          head_grid%write_metadata = .true.
---
>       !  We need to save the 3d data to compute a difference during the next loop.  Couple the
>       !  3d fields with total mu (mub + mu_2) and the stagger-specific map scale factor.
333,336c451,460
<          CALL calc_current_date ( grid%id , 0. )
<          head_grid%write_metadata = .true.
< !        CALL output_model_input ( id1, head_grid , config_flags , ierr )
< !        CALL close_dataset ( id1 , config_flags , "DATASET=INPUT" )
---
>       CALL couple ( grid%em_mu_2 , grid%em_mub , ubdy3dtemp1 , grid%em_u_2                 , 'u' , grid%msfu , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , vbdy3dtemp1 , grid%em_v_2                 , 'v' , grid%msfv , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , tbdy3dtemp1 , grid%em_t_2                 , 't' , grid%msft , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , pbdy3dtemp1 , grid%em_ph_2                , 'w' , grid%msft , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , qbdy3dtemp1 , grid%moist_2(:,:,:,P_QV) , 't' , grid%msft , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
338,339c462,466
<          !  Is this or is this not the last time time?  We can remove some unnecessary
<          !  stores if it is not.
---
>       DO j = jps , MIN(jde-1,jpe)
>          DO i = ips , MIN(ide-1,ipe)
>             mbdy2dtemp1(i,j) = grid%em_mu_2(i,j)
>          END DO
>       END DO
341c468,469
<          IF     ( loop .LT. time_loop_max ) THEN
---
>       !  There are 2 components to the lateral boundaries.  First, there is the starting
>       !  point of this time period - just the outer few rows and columns.
343,352c471,494
<             !  We need to save the 3d data to compute a difference during the next loop.  Couple the
<             !  3d fields with total mu (mub + mu_1) and the stagger-specific map scale factor.
<             !  We load up the boundary data again for use in the next loop.
<          
<             ubdy3dtemp1 = ubdy3dtemp2
<             vbdy3dtemp1 = vbdy3dtemp2
<             tbdy3dtemp1 = tbdy3dtemp2
<             pbdy3dtemp1 = pbdy3dtemp2
<             qbdy3dtemp1 = qbdy3dtemp2
<             mbdy2dtemp1 = mbdy2dtemp2
---
>       CALL stuff_bdy     ( ubdy3dtemp1 , grid%em_u_b     , 'U' , ijds , ijde , spec_bdy_width      , &
>                                                                  ids , ide , jds , jde , kds , kde , &
>                                                                  ims , ime , jms , jme , kms , kme , &
>                                                                  ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdy     ( vbdy3dtemp1 , grid%em_v_b     , 'V' , ijds , ijde , spec_bdy_width      , &
>                                                                  ids , ide , jds , jde , kds , kde , &
>                                                                  ims , ime , jms , jme , kms , kme , &
>                                                                  ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdy     ( tbdy3dtemp1 , grid%em_t_b     , 'T' , ijds , ijde , spec_bdy_width      , &
>                                                                  ids , ide , jds , jde , kds , kde , &
>                                                                  ims , ime , jms , jme , kms , kme , &
>                                                                  ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdy     ( pbdy3dtemp1 , grid%em_ph_b    , 'W' , ijds , ijde , spec_bdy_width      , &
>                                                                  ids , ide , jds , jde , kds , kde , &
>                                                                  ims , ime , jms , jme , kms , kme , &
>                                                                  ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdy     ( qbdy3dtemp1 , grid%em_rqv_b   , 'T' , ijds , ijde , spec_bdy_width      , &
>                                                                  ids , ide , jds , jde , kds , kde , &
>                                                                  ims , ime , jms , jme , kms , kme , &
>                                                                  ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdy2    ( mbdy2dtemp1 , grid%em_mu_b    , 'T' , ijds , ijde , spec_bdy_width      , &
>                                                                  ids , ide , jds , jde , kds , kde , &
>                                                                  ims , ime , jms , jme , kms , kme , &
>                                                                  ips , ipe , jps , jpe , kps , kpe )
354,368d495
<             !  There are 2 components to the lateral boundaries.  First, there is the starting
<             !  point of this time period - just the outer few rows and columns.
<       
<             CALL stuff_bdy     ( ubdy3dtemp1 , grid%em_u_b     , 'U' , ids , ide , jds , jde , kds , kde , &
<                                                                    ims , ime , jms , jme , kms , kme ) 
<             CALL stuff_bdy     ( vbdy3dtemp1 , grid%em_v_b     , 'V' , ids , ide , jds , jde , kds , kde , &
<                                                                    ims , ime , jms , jme , kms , kme ) 
<             CALL stuff_bdy     ( tbdy3dtemp1 , grid%em_t_b     , 'T' , ids , ide , jds , jde , kds , kde , &
<                                                                    ims , ime , jms , jme , kms , kme ) 
<             CALL stuff_bdy     ( pbdy3dtemp1 , grid%em_ph_b    , 'W' , ids , ide , jds , jde , kds , kde , &
<                                                                    ims , ime , jms , jme , kms , kme ) 
<             CALL stuff_bdy     ( qbdy3dtemp1 , grid%em_rqv_b   , 'T' , ids , ide , jds , jde , kds , kde , &
<                                                                    ims , ime , jms , jme , kms , kme ) 
<             CALL stuff_bdy2    ( mbdy2dtemp1 , grid%em_mu_b    , 'T' , ids , ide , jds , jde , kds , kde , &
<                                                                    ims , ime , jms , jme , kms , kme ) 
370c497
<          ELSE IF ( loop .EQ. time_loop_max ) THEN
---
>    ELSE IF ( loop .GT. 1 ) THEN
372c499
<             !  If this is the last time through here, we need to close the files.
---
>       !  Open the boundary file.
374d500
<             CALL close_dataset ( id , config_flags , "DATASET=BOUNDARY" )
376,377c502,512
<          END IF
< 
---
>       IF ( loop .eq. 2 ) THEN
>          CALL construct_filename1( bdyname , 'wrfbdy' , grid%id , 2 )
>          CALL open_w_dataset ( id, TRIM(bdyname) , grid , config_flags , output_boundary , "DATASET=BOUNDARY", ierr )
>          IF ( ierr .NE. 0 ) THEN
>                CALL wrf_error_fatal( 'real: error opening wrfbdy for writing' )
>          ENDIF
>          grid%write_metadata = .true.
>       ELSE
> ! what's this do?
>          grid%write_metadata = .false.
>          CALL ESMF_ClockAdvance( grid%domain_clock, rc=rc )
380d514
<       !  Process which time now?
382,384c516
<       current_date_char = current_date(1:19)
<       CALL geth_newdate ( next_date_char , current_date_char , interval_seconds ) 
<       start_date = next_date_char // '.0000'
---
>       !  Couple this time period's data with total mu, and save it in the *bdy3dtemp2 arrays.
386c518,527
<    END DO
---
>       CALL couple ( grid%em_mu_2 , grid%em_mub , ubdy3dtemp2 , grid%em_u_2                 , 'u' , grid%msfu , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , vbdy3dtemp2 , grid%em_v_2                 , 'v' , grid%msfv , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , tbdy3dtemp2 , grid%em_t_2                 , 't' , grid%msft , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , pbdy3dtemp2 , grid%em_ph_2                , 'w' , grid%msft , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
>       CALL couple ( grid%em_mu_2 , grid%em_mub , qbdy3dtemp2 , grid%moist_2(:,:,:,P_QV) , 't' , grid%msft , &
>                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
388c529,533
<    CALL       wrf_debug (   0 , 'wrf: SUCCESS COMPLETE REAL_EM INIT' )
---
>       DO j = jps , jpe
>          DO i = ips , ipe
>             mbdy2dtemp2(i,j) = grid%em_mu_2(i,j)
>          END DO
>       END DO
390,392c535,537
< #ifdef DM_PARALLEL
<     CALL wrf_dm_shutdown
< #endif
---
>       !  During all of the loops after the first loop, we first compute the boundary
>       !  tendencies with the current data values (*bdy3dtemp2 arrays) and the previously 
>       !  saved information stored in the *bdy3dtemp1 arrays.
394c539,673
< END PROGRAM real_data
---
>       CALL stuff_bdytend ( ubdy3dtemp2 , ubdy3dtemp1 , REAL(interval_seconds) , grid%em_u_bt  , 'U' , &
>                                                             ijds , ijde , spec_bdy_width      , &
>                                                             ids , ide , jds , jde , kds , kde , &
>                                                             ims , ime , jms , jme , kms , kme , &
>                                                             ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdytend ( vbdy3dtemp2 , vbdy3dtemp1 , REAL(interval_seconds) , grid%em_v_bt  , 'V' , &
>                                                             ijds , ijde , spec_bdy_width      , &
>                                                             ids , ide , jds , jde , kds , kde , &
>                                                             ims , ime , jms , jme , kms , kme , &
>                                                             ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdytend ( tbdy3dtemp2 , tbdy3dtemp1 , REAL(interval_seconds) , grid%em_t_bt  , 'T' , &
>                                                             ijds , ijde , spec_bdy_width      , &
>                                                             ids , ide , jds , jde , kds , kde , &
>                                                             ims , ime , jms , jme , kms , kme , &
>                                                             ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdytend ( pbdy3dtemp2 , pbdy3dtemp1 , REAL(interval_seconds) , grid%em_ph_bt  , 'W' , &
>                                                             ijds , ijde , spec_bdy_width      , &
>                                                             ids , ide , jds , jde , kds , kde , &
>                                                             ims , ime , jms , jme , kms , kme , &
>                                                             ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdytend ( qbdy3dtemp2 , qbdy3dtemp1 , REAL(interval_seconds) , grid%em_rqv_bt , 'T' , &
>                                                             ijds , ijde , spec_bdy_width      , &
>                                                             ids , ide , jds , jde , kds , kde , &
>                                                             ims , ime , jms , jme , kms , kme , &
>                                                             ips , ipe , jps , jpe , kps , kpe )
>       CALL stuff_bdytend2( mbdy2dtemp2 , mbdy2dtemp1 , REAL(interval_seconds) , grid%em_mu_bt  , 'T' , &
>                                                             ijds , ijde , spec_bdy_width      , &
>                                                             ids , ide , jds , jde , kds , kde , &
>                                                             ims , ime , jms , jme , kms , kme , &
>                                                             ips , ipe , jps , jpe , kps , kpe )
> 
>       !  Both pieces of the boundary data are now available to be written (initial time and tendency).
>       !  This looks ugly, these date shifting things.  What's it for?  We want the "Times" variable
>       !  in the lateral BDY file to have the valid times of when the initial fields are written.
>       !  That's what the loop-2 thingy is for with the start date.  We increment the start_date so
>       !  that the starting time in the attributes is the second time period.  Why you may ask.  I
>       !  agree, why indeed.
> 
>       CALL atotime( current_date(1:19), grid%current_time )
>       CALL ESMF_ClockSetCurrTime(grid%domain_clock, grid%current_time, rc)
> 
>       temp24= current_date
>       temp24b=start_date
>       start_date = current_date
>       CALL geth_newdate ( temp19 , temp24b(1:19) , (loop-2) * model_config_rec%interval_seconds )
>       current_date = temp19 //  '.0000'
>       CALL atotime( current_date(1:19), grid%current_time )
>       CALL ESMF_ClockSetCurrTime(grid%domain_clock, grid%current_time, rc)
> print *,'LBC valid between these times ',current_date, ' ',start_date
>       CALL output_boundary ( id, grid , config_flags , ierr )
>       current_date = temp24
>       start_date = temp24b
>       CALL atotime( current_date(1:19), grid%current_time )
>       CALL ESMF_ClockSetCurrTime(grid%domain_clock, grid%current_time, rc)
> 
>       !  OK, for all of the loops, we output the initialzation data, which would allow us to
>       !  start the model at any of the available analysis time periods.
> 
> !     WRITE ( loop_char , FMT = '(I4.4)' ) loop
> !     CALL open_w_dataset ( id1, 'wrfinput'//loop_char , grid , config_flags , output_model_input , "DATASET=INPUT", ierr )
> !     IF ( ierr .NE. 0 ) THEN
> !       CALL wrf_error_fatal( 'real: error opening wrfinput'//loop_char//' for writing' )
> !     ENDIF
> !     grid%write_metadata = .true.
> 
> !     CALL calc_current_date ( grid%id , 0. )
> !     CALL output_model_input ( id1, grid , config_flags , ierr )
> !     CALL close_dataset ( id1 , config_flags , "DATASET=INPUT" )
> 
>       !  Is this or is this not the last time time?  We can remove some unnecessary
>       !  stores if it is not.
> 
>       IF     ( loop .LT. time_loop_max ) THEN
> 
>          !  We need to save the 3d data to compute a difference during the next loop.  Couple the
>          !  3d fields with total mu (mub + mu_2) and the stagger-specific map scale factor.
>          !  We load up the boundary data again for use in the next loop.
> 
>          DO j = jps , jpe
>             DO k = kps , kpe
>                DO i = ips , ipe
>                   ubdy3dtemp1(i,k,j) = ubdy3dtemp2(i,k,j)
>                   vbdy3dtemp1(i,k,j) = vbdy3dtemp2(i,k,j)
>                   tbdy3dtemp1(i,k,j) = tbdy3dtemp2(i,k,j)
>                   pbdy3dtemp1(i,k,j) = pbdy3dtemp2(i,k,j)
>                   qbdy3dtemp1(i,k,j) = qbdy3dtemp2(i,k,j)
>                END DO
>             END DO
>          END DO
> 
>          DO j = jps , jpe
>             DO i = ips , ipe
>                mbdy2dtemp1(i,  j) = mbdy2dtemp2(i,  j)
>             END DO
>          END DO
> 
>          !  There are 2 components to the lateral boundaries.  First, there is the starting
>          !  point of this time period - just the outer few rows and columns.
> 
>          CALL stuff_bdy     ( ubdy3dtemp1 , grid%em_u_b     , 'U' , ijds , ijde , spec_bdy_width      , &
>                                                                     ids , ide , jds , jde , kds , kde , &
>                                                                     ims , ime , jms , jme , kms , kme , &
>                                                                     ips , ipe , jps , jpe , kps , kpe )
>          CALL stuff_bdy     ( vbdy3dtemp1 , grid%em_v_b     , 'V' , ijds , ijde , spec_bdy_width      , &
>                                                                     ids , ide , jds , jde , kds , kde , &
>                                                                     ims , ime , jms , jme , kms , kme , &
>                                                                     ips , ipe , jps , jpe , kps , kpe )
>          CALL stuff_bdy     ( tbdy3dtemp1 , grid%em_t_b     , 'T' , ijds , ijde , spec_bdy_width      , &
>                                                                     ids , ide , jds , jde , kds , kde , &
>                                                                     ims , ime , jms , jme , kms , kme , &
>                                                                     ips , ipe , jps , jpe , kps , kpe )
>          CALL stuff_bdy     ( pbdy3dtemp1 , grid%em_ph_b    , 'W' , ijds , ijde , spec_bdy_width      , &
>                                                                     ids , ide , jds , jde , kds , kde , &
>                                                                     ims , ime , jms , jme , kms , kme , &
>                                                                     ips , ipe , jps , jpe , kps , kpe )
>          CALL stuff_bdy     ( qbdy3dtemp1 , grid%em_rqv_b   , 'T' , ijds , ijde , spec_bdy_width      , &
>                                                                     ids , ide , jds , jde , kds , kde , &
>                                                                     ims , ime , jms , jme , kms , kme , &
>                                                                     ips , ipe , jps , jpe , kps , kpe )
>          CALL stuff_bdy2    ( mbdy2dtemp1 , grid%em_mu_b    , 'T' , ijds , ijde , spec_bdy_width      , &
>                                                                     ids , ide , jds , jde , kds , kde , &
>                                                                     ims , ime , jms , jme , kms , kme , &
>                                                                     ips , ipe , jps , jpe , kps , kpe )
> 
>       ELSE IF ( loop .EQ. time_loop_max ) THEN
> 
>          !  If this is the last time through here, we need to close the files.
> 
>          CALL close_dataset ( id , config_flags , "DATASET=BOUNDARY" )
> 
>       END IF
> 
>    END IF
> 
> END SUBROUTINE assemble_output