tem_unify_surfaceData Subroutine

  subroutine tem_unify_surfaceData( me, all_pointCoords )
    ! ---------------------------------------------------------------------------
    !> datatype to store the surface information
    type( tem_surfData_type ), intent(inout) :: me
    !> tmp point coordinates to be unified and stored in me
    real(kind=rk), intent(in) :: all_pointCoords(:,:)
    ! ---------------------------------------------------------------------------
    ! temporary map from the all_pointCoords to the unique ones
    integer :: map2unique(size(all_pointCoords,2))
    ! array of unique point coordinates in physical coordinates
    ! size: 3, nUniquePoints_total
    type( grw_real2darray_type ) :: pointCoords
    ! temporary dynamic array of unique treeIDs (for unifying the points)
    type( dyn_longArray_type ) :: uniqueTreeIDs
    ! growing array of integers counting the number of identical points
    type( grw_intArray_type ) :: grw_counter
    ! temporary coordinates of an element (iX,iY,iZ, level)
    integer :: tmp_coord(4)
    ! temporary treeID
    integer(kind=long_k) :: treeID
    ! position of the element in the dynamic array of unique treeIDs, in the
    ! growing array of PointCoords and in the growing array of counters
    integer :: pos
    ! was this treeID appended to the unique list
    logical :: wasAdded
    ! local tree with no treeIDs needed for unification of the points
    type( treelmesh_type ) :: loc_tree
    ! minimum x, y, z coordinate
    real(kind=rk) :: minX, minY, minZ
    ! maximum x, y, z coordinate
    real(kind=rk) :: maxX, maxY, maxZ
    ! counters
    integer :: iCoord
    integer :: iTria
    ! ---------------------------------------------------------------------------

    write(logUnit(10),*) " Unifying the coordinates ..."

    ! initialize the growing and dynamic arrays
    call init( me = pointCoords, width = 3 )

    ! get the min and max of the coordinates
    minX = minval(all_pointCoords(1,:))
    minY = minval(all_pointCoords(2,:))
    minZ = minval(all_pointCoords(3,:))

    maxX = maxval(all_pointCoords(1,:))
    maxY = maxval(all_pointCoords(2,:))
    maxZ = maxval(all_pointCoords(3,:))

    ! initialize the local tree with an origin and a bc_length
    loc_tree%global%Origin = (/minX,minY,minZ/)
    loc_tree%global%BoundingCubeLength = max( maxX - minX,                     &
      &                                       maxY - minY,                     &
      &                                       maxZ - minZ )

    ! get the treeIDs on the highest level to unify the coordinates and map
    ! the coordinates on these unique treeIDs
    do iCoord = 1, size( all_pointCoords, 2)
      ! get the coordinates
      tmp_coord = tem_CoordOfReal( loc_tree, all_pointCoords(:,iCoord ) )
      treeID = tem_IdOfCoord( tmp_coord )
      call append( me       = uniqueTreeIDs,                                   &
        &          val      = treeID,                                          &
        &          pos      = pos,                                             &
        &          wasAdded = wasAdded )
      map2unique(iCoord) = pos
      ! if the value was added the first time ...
      if( wasAdded )then
        ! ... append the coordinate
        call append( me  = pointCoords,                                        &
          &          val = all_pointCoords(:,iCoord ),                         &
          &          pos = pos )
        ! ... set the counter to 1 for this element
        call placeAt( me  = grw_counter, &
          &           val = 1,           &
          &           pos = pos          )
      else ! if the value was existing before ...
        ! ...read the coordinate and multiply it with the counter (number of
        ! times this element was added before) add the new coordinate average
        ! it by counter + 1 and write it back
        pointCoords%val( :, pos ) = ( pointCoords%val( :, pos )                &
          &     * grw_counter%val( pos ) + all_pointCoords( :, iCoord ))       &
          &     / ( grw_counter%val( pos ) + 1 )
        ! ... and increase the counter by one
        grw_counter%val( pos ) = grw_counter%val( pos ) + 1
      end if
    end do

    ! initialize the number of unique points and allocate the pointCoords
    ! array
    me%nUniquePoints_total = pointCoords%nVals
    allocate( me%pointCoords( 3*me%nUniquePoints_total ))
    do iCoord = 1, me%nUniquePoints_total
      me%PointCoords((iCoord-1)*3+1:(iCoord-1)*3+3) = pointCoords%val(:,iCoord)
    end do

    ! loop over the triangles and triangle vertices and store the right
    ! coordinate (barycenter of the parent element on the highest level)
    ! and update the position in the trias connectivity array
    do iTria = 1, me%nTrias
      do iCoord = 1, 3
        ! update the position in the triangle connectivity array
        me%trias( iCoord, iTria ) = map2unique( me%trias( iCoord, iTria ))
      end do
    end do

    call destroy( uniqueTreeIDs )
    call destroy( grw_counter )
    call destroy( pointCoords )

    write(logUnit(1),*) " nTriangles: ", me%nTrias
    write(logUnit(1),*) " nPoints: ", me%nUniquePoints_total

  end subroutine tem_unify_surfaceData