如何在R中高效处理RasterBrick实现日均温年度分布离散化？

我手里有一个存储德国15年历史日均温的大型RasterBrick栅格数据，数据基本信息如下：

> Temperature_rasterData 
class      : RasterBrick 
dimensions : 31, 37, 1147, 5479  (nrow, ncol, ncell, nlayers)
resolution : 0.25, 0.25  (x, y)
extent     : 5.75, 15, 47.25, 55  (xmin, xmax, ymin, ymax)
coord. ref.: +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source: in memory
names      : X1980.01.01, X1980.01.02, X1980.01.03, X1980.01.04, X1980.01.05, X1980.01.06, X1980.01.07, X1980.01.08, X1980.01.09, X1980.01.10, X1980.01.11, X1980.01.12, X1980.01.13, X1980.01.14, X1980.01.15, ... 
min values : -9.24, -11.32, -12.05, -14.12, -7.91, -6.35, -6.74, -7.77, -9.79, -10.17, -12.20, -14.90, -15.68, -15.61, -15.22, ... 
max values : 2.19, 0.68, 0.30, 2.91, 5.25, 5.03, 4.33, 3.40, 1.52, 0.33, -1.10, -1.61, -3.55, -0.12, 0.19, ... 

我的核心需求是：把日均温的年度分布离散化为10个固定温度区间（参考相关文献的温度区间划分方法），具体步骤需要先获取所有栅格图层的温度最大/最小值来确定整体温度范围，再对每个网格的年度日均温分布做区间划分。

因为原始数据体量太大没法直接复现，我准备了小型测试栅格数据，处理这个测试数据的R脚本如下：

temp_raster <- raster::stack('~/tg_day_2017_grid_ensmean.nc')
data(wrld_simpl)
Germany <- wrld_simpl[wrld_simpl@data$NAME == "Germany",]
deu_ext <- extent(Germany)
Deu_crop <- crop(temp_raster ,deu_ext)

之前我试过用下面的代码获取温度范围，但总觉得这个方法不够简洁高效，处理大型数据的时候应该会很慢：

nms <- names(Deu_crop)
yrs <- unique(sub('X(\\d+).+','\\1',nms))
getRange <- lapply(yrs,function(x) { 
  range(Deu_crop[[grep(x,nms)]],na.rm=TRUE) 
})

现在我卡壳的地方主要有两个：

1. 不知道怎么在大型RasterBrick里高效实现数据离散化，尤其是针对这种包含上千个日均温观测图层的栅格数据；
2. 希望能按德国NUTS区域（我已经准备好对应的Shapefile）对每年的日均温分布进行离散化统计，并且绘制类似下图的年度柱状图：

恳请各位大佬分享一下在R中实现这些需求的简便方法，万分感谢！

内容的提问来源于stack exchange，提问作者Andy.Jian