| 1 | subroutine MPP_READ_2DDECOMP_1D_( unit, field, domain, data, tindex ) |
|---|
| 2 | integer, intent(in) :: unit |
|---|
| 3 | type(fieldtype), intent(in) :: field |
|---|
| 4 | type(domain2D), intent(in) :: domain |
|---|
| 5 | MPP_TYPE_, intent(inout) :: data(:) |
|---|
| 6 | integer, intent(in), optional :: tindex |
|---|
| 7 | integer::il_xbegin,il_xend,il_ybegin,il_yend |
|---|
| 8 | MPP_TYPE_,allocatable :: data3D(:,:,:) |
|---|
| 9 | |
|---|
| 10 | call mpp_get_compute_domain(domain & |
|---|
| 11 | ,xbegin=il_xbegin,xend=il_xend & |
|---|
| 12 | ,ybegin=il_ybegin,yend=il_yend) |
|---|
| 13 | |
|---|
| 14 | allocate(data3D(1:il_xend-il_xbegin+1,1:il_yend-il_ybegin+1,1:1)) |
|---|
| 15 | data3D = RESHAPE( data, SHAPE(data3D) ) |
|---|
| 16 | call mpp_read( unit, field, domain, data3D, tindex ) |
|---|
| 17 | data = RESHAPE( data3D, SHAPE(data) ) |
|---|
| 18 | deallocate(data3D) |
|---|
| 19 | |
|---|
| 20 | return |
|---|
| 21 | end subroutine MPP_READ_2DDECOMP_1D_ |
|---|
| 22 | subroutine MPP_READ_2DDECOMP_2D_( unit, field, domain, data, tindex ) |
|---|
| 23 | integer, intent(in) :: unit |
|---|
| 24 | type(fieldtype), intent(in) :: field |
|---|
| 25 | type(domain2D), intent(in) :: domain |
|---|
| 26 | MPP_TYPE_, intent(inout) :: data(:,:) |
|---|
| 27 | integer, intent(in), optional :: tindex |
|---|
| 28 | MPP_TYPE_ :: data3D(size(data,1),size(data,2),1) |
|---|
| 29 | #ifdef use_CRI_pointers |
|---|
| 30 | pointer( ptr, data3D ) |
|---|
| 31 | ptr = LOC(data) |
|---|
| 32 | call mpp_read( unit, field, domain, data3D, tindex ) |
|---|
| 33 | #else |
|---|
| 34 | data3D = RESHAPE( data, SHAPE(data3D) ) |
|---|
| 35 | call mpp_read( unit, field, domain, data3D, tindex ) |
|---|
| 36 | data = RESHAPE( data3D, SHAPE(data) ) |
|---|
| 37 | #endif |
|---|
| 38 | return |
|---|
| 39 | end subroutine MPP_READ_2DDECOMP_2D_ |
|---|
| 40 | |
|---|
| 41 | subroutine MPP_READ_2DDECOMP_3D_( unit, field, domain, data, tindex ) |
|---|
| 42 | !mpp_read reads <data> which has the domain decomposition <domain> |
|---|
| 43 | integer, intent(in) :: unit |
|---|
| 44 | type(fieldtype), intent(in) :: field |
|---|
| 45 | type(domain2D), intent(in) :: domain |
|---|
| 46 | MPP_TYPE_, intent(inout) :: data(:,:,:) |
|---|
| 47 | integer, intent(in), optional :: tindex |
|---|
| 48 | MPP_TYPE_, allocatable :: cdata(:,:,:) |
|---|
| 49 | MPP_TYPE_, allocatable :: gdata(:) |
|---|
| 50 | integer :: len, lenx,leny,lenz,i,j,k,n |
|---|
| 51 | !NEW: data may be on compute OR data domain |
|---|
| 52 | logical :: data_has_halos, halos_are_global, x_is_global, y_is_global |
|---|
| 53 | integer :: is, ie, js, je, isd, ied, jsd, jed, isg, ieg, jsg, jeg, ioff, joff |
|---|
| 54 | |
|---|
| 55 | if (.NOT. present(tindex) .AND. mpp_file(unit)%time_level .ne. -1) & |
|---|
| 56 | call mpp_error(FATAL, 'MPP_READ: need to specify a time level for data with time axis') |
|---|
| 57 | |
|---|
| 58 | if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_READ: must first call mpp_io_init.' ) |
|---|
| 59 | if( .NOT.mpp_file(unit)%opened )call mpp_error( FATAL, 'MPP_READ: invalid unit number.' ) |
|---|
| 60 | |
|---|
| 61 | call mpp_get_compute_domain( domain, is, ie, js, je ) |
|---|
| 62 | call mpp_get_data_domain ( domain, isd, ied, jsd, jed, x_is_global=x_is_global, y_is_global=y_is_global ) |
|---|
| 63 | call mpp_get_global_domain ( domain, isg, ieg, jsg, jeg ) |
|---|
| 64 | if( size(data,1).EQ.ie-is+1 .AND. size(data,2).EQ.je-js+1 )then |
|---|
| 65 | data_has_halos = .FALSE. |
|---|
| 66 | else if( size(data,1).EQ.ied-isd+1 .AND. size(data,2).EQ.jed-jsd+1 )then |
|---|
| 67 | data_has_halos = .TRUE. |
|---|
| 68 | else |
|---|
| 69 | call mpp_error( FATAL, 'MPP_READ: data must be either on compute domain or data domain.' ) |
|---|
| 70 | end if |
|---|
| 71 | halos_are_global = x_is_global .AND. y_is_global |
|---|
| 72 | if( npes.GT.1 .AND. mpp_file(unit)%threading.EQ.MPP_SINGLE )then |
|---|
| 73 | if( halos_are_global )then !you can read directly into data array |
|---|
| 74 | if( pe.EQ.0 )call read_record( unit, field, size(data), data, tindex ) |
|---|
| 75 | else |
|---|
| 76 | lenz=size(data,3) |
|---|
| 77 | len=lenx*leny*lenz |
|---|
| 78 | allocate(gdata(len)) |
|---|
| 79 | ! read field on pe 0 and pass to all pes |
|---|
| 80 | if( pe.EQ.0 ) call read_record( unit, field, len, gdata, tindex ) |
|---|
| 81 | ! broadcasting global array, this can be expensive! |
|---|
| 82 | call mpp_transmit( gdata, len, ALL_PES, gdata, len, 0 ) |
|---|
| 83 | ioff = is; joff = js |
|---|
| 84 | if( data_has_halos )then |
|---|
| 85 | ioff = isd; joff = jsd |
|---|
| 86 | end if |
|---|
| 87 | do k=1,size(data,3) |
|---|
| 88 | do j=js,je |
|---|
| 89 | do i=is,ie |
|---|
| 90 | n=(i-isg+1) + (j-jsg)*lenx + (k-1)*lenx*leny |
|---|
| 91 | data(i-ioff+1,j-joff+1,k)=gdata(n) |
|---|
| 92 | enddo |
|---|
| 93 | enddo |
|---|
| 94 | enddo |
|---|
| 95 | deallocate(gdata) |
|---|
| 96 | end if |
|---|
| 97 | else if( data_has_halos )then |
|---|
| 98 | ! for uniprocessor or multithreaded read |
|---|
| 99 | ! read compute domain as contiguous data |
|---|
| 100 | |
|---|
| 101 | allocate( cdata(is:ie,js:je,size(data,3)) ) |
|---|
| 102 | call read_record(unit,field,size(cdata),cdata,tindex,domain) |
|---|
| 103 | |
|---|
| 104 | data(is-isd+1:ie-isd+1,js-jsd+1:je-jsd+1,:) = cdata(:,:,:) |
|---|
| 105 | deallocate(cdata) |
|---|
| 106 | else |
|---|
| 107 | call read_record(unit,field,size(data),data,tindex,domain) |
|---|
| 108 | end if |
|---|
| 109 | |
|---|
| 110 | return |
|---|
| 111 | end subroutine MPP_READ_2DDECOMP_3D_ |
|---|
