1Glacier-level calculation.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "32d01cad-de27-47a3-9ae0-a420bb1d63ca",
   "metadata": {},
   "outputs": [],
   "source": [
    "import geopandas as gpd\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import os\n",
    "import xarray as xr\n",
    "import rioxarray as rioxr\n",
    "import numpy as np\n",
    "import xdem\n",
    "import geoutils as gu\n",
    "import warnings\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "fb52a7d9-2f1a-48d3-94c2-f710d7e0e3f0",
   "metadata": {},
   "outputs": [],
   "source": [
    "pth='D:\\\\Key_works\\\\20231212dh_calculate'\n",
    "pth_gla='D:\\\\Key_works\\\\20231212dh_calculate\\\\gla\\\\CGI1'\n",
    "pth_srtm='D:\\\\Key_works\\\\20231212dh_calculate\\\\srtm'\n",
    "pth_srtm_fill='D:\\\\Key_works\\\\20231212dh_calculate\\\\srtm_fill'\n",
    "pth_dh_xnhl=pth+'/xnhl_test/dh_pene_all'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "76d22ec9-d239-4cc3-80d8-49829979b691",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "261"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dxt_ix=pth+'/xnhl_test/shp/dxt_group_xn.shp'\n",
    "df_dxt=gpd.read_file(dxt_ix)\n",
    "####选取研究区内的dh\n",
    "group_lst=[]\n",
    "for i in df_dxt.reg_group.values:\n",
    "    if i is not None:\n",
    "        ist=i.split('+')\n",
    "        for j in ist:\n",
    "            group_lst.append(j)\n",
    "group_lst=list(set(group_lst))# 去除重复元素 保证标识的唯一性\n",
    "len(group_lst)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "12e03bec-2481-4e22-a54e-ec9d6f900c14",
   "metadata": {},
   "outputs": [],
   "source": [
    "def nc_out_stp1(gla_dh,gla_dh_fill,gla_h,gpi,n_num,error,outnm_nc):\n",
    "    ########输出topo dh 的netcdf版本#############################\n",
    "    mytrans=gla_h.coords['spatial_ref'].GeoTransform\n",
    "    mytrans=[float(i) for i in mytrans.split(' ')]\n",
    "    myattrs={'transform': (mytrans[1], mytrans[2], mytrans[0], mytrans[4], mytrans[5], mytrans[3]),\n",
    "             'crs': '+proj=aea +lat_0=0 +lon_0=105 +lat_1=25 +lat_2=47 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs=True',\n",
    "             'res': (30.0, 30.0),'nodatavals': 'nan','scales': (1.0,),'offsets': (0.0,),\n",
    "             'topo process group': '%s'%gpi,\n",
    "             'effective samples':round(float(n_num),4),\n",
    "             'dh error':round(float(error),4),}\n",
    "    ds_nc_out=xr.Dataset(data_vars={'topo_dh':(['y','x'],gla_dh.values[0]),\n",
    "                                    'elevation':(['y','x'],gla_h.values[0]),\n",
    "                                    'topo_dh_fill':(['y','x'],gla_dh_fill.values[0])\n",
    "                                   \n",
    "                                   },\n",
    "                         coords={'y':gla_dh.y.values,'x':gla_dh.x.values},attrs=myattrs)\n",
    "    ds_nc_out['topo_dh'].attrs={'long_name':'Elevation difference from Topo DEM and SRTM DEM (viod version, v1)','units':'m'}\n",
    "    ds_nc_out['topo_dh_fill'].attrs={'long_name':'Filled elevation difference from Topo DEM and SRTM DEM','units':'m','Fill method':'Normalized regional hypsometric interpolation'}\n",
    "    ds_nc_out['elevation'].attrs={'long_name':'elevation from viod filled SRTM DEM','units':'m a.s.l'}\n",
    "    ds_nc_out.to_netcdf(outnm_nc)\n",
    "\n",
    "###处理DH中异常值  flag=1表示处理\n",
    "##flag=0表示小冰川不处理\n",
    "def outlier_rm_ip(gla_dh0,gla_h):\n",
    "    gla_dh=gla_dh0.where((gla_dh0>-150)&(gla_dh0<150))\n",
    "    ####小冰川不处理，处理之后空值插值拟合有问题\n",
    "    if gla_dh.count().values<10:\n",
    "        return gla_dh,0,gla_dh.count().values/gla_h.count().values\n",
    "    #############################################################\n",
    "    gla_ma=gla_h.median().values  #物质平衡线高度，RGI中有些是错误的，因此，这里直接计算\n",
    "    mask=gla_h.values*0+1  #冰川区mask  \n",
    "    mask_acc=(gla_h.where(gla_h>=gla_ma)).values*0+1###积累区mask\n",
    "    mask_abla=(gla_h.where(gla_h<gla_ma)).values*0+1 ##消融区mask\n",
    "    pos1=np.where(gla_h.values>=gla_ma) ###积累区位置\n",
    "    pos2=np.where(gla_h.values<gla_ma) ###消融区位置\n",
    "    #############################################################\n",
    "    dh_abla=gla_dh*mask_abla\n",
    "    ##################outlier remove  ref. Zhou 2018 RSE\n",
    "    w=(gla_h.max()-gla_h)/(gla_h.max()-gla_h.min())\n",
    "    A=abs(dh_abla).max().values*2\n",
    "    h_delta=A*w**2\n",
    "    ggs1=gla_dh*1\n",
    "    msst=(abs(gla_dh)-abs(h_delta)) \n",
    "    #h_delta.plot()\n",
    "    #print(h_delta)\n",
    "    pps1=np.where(msst.values>0)##剔除区\n",
    "    gla_dh.values[pps1]=np.nan\n",
    "    ###################################################\n",
    "    dh_acc=gla_dh*mask_acc\n",
    "    dh_abla=gla_dh*mask_abla\n",
    "\n",
    "    dh_acc_med=dh_acc.median().values\n",
    "    dh_abla_med=dh_abla.median().values\n",
    "           \n",
    "    return gla_dh,1,gla_dh.count().values/gla_h.count().values\n",
    "\n",
    "\n",
    "'''generate strip positions and interval center by interval'''\n",
    "def ExtractByStrip(df,strip_sign,outsign,interval):\n",
    "    #val_min=df[strip_sign].min()\n",
    "    #val_max=df[strip_sign].max()\n",
    "    val_min=1000\n",
    "    val_max=8500\n",
    "    #print(val_min,val_max)\n",
    "    #n=int((val_max-val_min)/interval)\n",
    "    df[outsign]=-1\n",
    "    x=val_min\n",
    "    while x<val_max:\n",
    "        y=x+interval\n",
    "        if y>=val_max:\n",
    "            y=val_max+1\n",
    "        df.loc[((df[strip_sign]>=x)&(df[strip_sign]<y)),outsign] =(x+y)/2.0\n",
    "        x=y\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "cd6ff975-653e-4652-9e30-57a4831f9553",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "g0231 done!\n",
      "g0111 done!\n",
      "g0200 done!\n",
      "g0265 done!\n",
      "g0337 done!\n",
      "g0365 done!\n",
      "g0185 done!\n",
      "g0123 done!\n",
      "g0273 done!\n",
      "g0199 done!\n",
      "g0197 done!\n",
      "g0158 done!\n",
      "g0150 done!\n",
      "g0175 done!\n",
      "g0027 done!\n",
      "g0334 done!\n",
      "g0114 done!\n",
      "g0100 done!\n",
      "g0173 done!\n",
      "g0217 done!\n",
      "g0232 done!\n",
      "g0151 done!\n",
      "g0033 done!\n",
      "g0209 done!\n",
      "g0349 done!\n",
      "g0176 done!\n",
      "g0017 done!\n",
      "g0061 done!\n",
      "g0041 done!\n",
      "g0004 done!\n",
      "g0142 done!\n",
      "g0051 done!\n",
      "g0312 done!\n",
      "g0240 done!\n",
      "g0152 done!\n",
      "g0331 done!\n",
      "g0011 done!\n",
      "g0424 done!\n",
      "g0126 done!\n",
      "g0120 done!\n",
      "g0174 done!\n",
      "g0347 done!\n",
      "g0196 done!\n",
      "g0044 done!\n",
      "g0178 done!\n",
      "g0329 done!\n",
      "g0108 done!\n",
      "g0235 done!\n",
      "g0224 done!\n",
      "g0201 done!\n",
      "g0256 done!\n",
      "g0190 done!\n",
      "g0162 done!\n",
      "g0202 done!\n",
      "g0259 done!\n",
      "g0211 done!\n",
      "g0030 done!\n",
      "g0062 done!\n",
      "g0106 done!\n",
      "g0253 done!\n",
      "g0335 done!\n",
      "g0022 done!\n",
      "g0155 done!\n",
      "g0109 done!\n",
      "g0043 done!\n",
      "g0169 done!\n",
      "g0166 done!\n",
      "g0052 done!\n",
      "g0180 done!\n",
      "g0286 done!\n",
      "g0188 done!\n",
      "g0414 done!\n",
      "g0278 done!\n",
      "g0270 done!\n",
      "g0315 done!\n",
      "g0136 done!\n",
      "g0204 done!\n",
      "g0317 done!\n",
      "g0034 done!\n",
      "g0262 done!\n",
      "g0219 done!\n",
      "g0345 done!\n",
      "g0213 done!\n",
      "g0157 done!\n",
      "g0145 done!\n",
      "g0326 done!\n",
      "g0220 done!\n",
      "g0124 done!\n",
      "g0177 done!\n",
      "g0333 done!\n",
      "g0026 done!\n",
      "g0121 done!\n",
      "g0229 done!\n",
      "g0206 done!\n",
      "g0218 done!\n",
      "g0128 done!\n",
      "g0319 done!\n",
      "g0130 done!\n",
      "g0322 done!\n",
      "g0007 done!\n",
      "g0104 done!\n",
      "g0098 done!\n",
      "g0116 done!\n",
      "g0129 done!\n",
      "g0343 done!\n",
      "g0141 done!\n",
      "g0183 done!\n",
      "g0346 done!\n",
      "g0258 done!\n",
      "g0137 done!\n",
      "g0113 done!\n",
      "g0321 done!\n",
      "g0182 done!\n",
      "g0005 done!\n",
      "g0159 done!\n",
      "g0257 done!\n",
      "g0226 done!\n",
      "g0274 done!\n",
      "g0138 done!\n",
      "g0035 done!\n",
      "g0214 done!\n",
      "g0281 done!\n",
      "g0189 done!\n",
      "g0163 done!\n",
      "g0046 done!\n",
      "g0283 done!\n",
      "g0110 done!\n",
      "g0016 done!\n",
      "g0195 done!\n",
      "g0203 done!\n",
      "g0288 done!\n",
      "g0342 done!\n",
      "g0002 done!\n",
      "g0134 done!\n",
      "g0023 done!\n",
      "g0193 done!\n",
      "g0207 done!\n",
      "g0267 done!\n",
      "g0097 done!\n",
      "g0029 done!\n",
      "g0228 done!\n",
      "g0179 done!\n",
      "g0413 done!\n",
      "g0341 done!\n",
      "g0271 done!\n",
      "g0316 done!\n",
      "g0132 done!\n",
      "g0127 done!\n",
      "g0119 done!\n",
      "g0149 done!\n",
      "g0167 done!\n",
      "g0006 done!\n",
      "g0059 done!\n",
      "g0101 done!\n",
      "g0105 done!\n",
      "g0282 done!\n",
      "g0147 done!\n",
      "g0172 done!\n",
      "g0324 done!\n",
      "g0058 done!\n",
      "g0140 done!\n",
      "g0252 done!\n",
      "g0001 done!\n",
      "g0031 done!\n",
      "g0323 done!\n",
      "g0212 done!\n",
      "g0170 done!\n",
      "g0118 done!\n",
      "g0148 done!\n",
      "g0165 done!\n",
      "g0423 done!\n",
      "g0239 done!\n",
      "g0171 done!\n",
      "g0332 done!\n",
      "g0330 done!\n",
      "g0208 done!\n",
      "g0025 done!\n",
      "g0336 done!\n",
      "g0210 done!\n",
      "g0215 done!\n",
      "g0225 done!\n",
      "g0103 done!\n",
      "g0221 done!\n",
      "g0205 done!\n",
      "g0125 done!\n",
      "g0313 done!\n",
      "g0366 done!\n",
      "g0416 done!\n",
      "g0272 done!\n",
      "g0261 done!\n",
      "g0249 done!\n",
      "g0181 done!\n",
      "g0133 done!\n",
      "g0143 done!\n",
      "g0216 done!\n",
      "g0234 done!\n",
      "g0168 done!\n",
      "g0039 done!\n",
      "g0233 done!\n",
      "g0194 done!\n",
      "g0340 done!\n",
      "g0328 done!\n",
      "g0339 done!\n",
      "g0250 done!\n",
      "g0015 done!\n",
      "g0238 done!\n",
      "g0131 done!\n",
      "g0112 done!\n",
      "g0045 done!\n",
      "g0018 done!\n",
      "g0008 done!\n",
      "g0107 done!\n",
      "g0038 done!\n",
      "g0161 done!\n",
      "g0020 done!\n",
      "g0117 done!\n",
      "g0191 done!\n",
      "g0338 done!\n",
      "g0060 done!\n",
      "g0236 done!\n",
      "g0032 done!\n",
      "g0102 done!\n",
      "g0320 done!\n",
      "g0318 done!\n",
      "g0192 done!\n",
      "g0348 done!\n",
      "g0309 done!\n",
      "g0156 done!\n",
      "g0251 done!\n",
      "g0146 done!\n",
      "g0040 done!\n",
      "g0344 done!\n",
      "g0230 done!\n",
      "g0237 done!\n",
      "g0024 done!\n",
      "g0184 done!\n",
      "g0285 done!\n",
      "g0144 done!\n",
      "g0139 done!\n",
      "g0160 done!\n",
      "g0222 done!\n",
      "g0255 done!\n",
      "g0135 done!\n",
      "g0198 done!\n",
      "g0311 done!\n",
      "g0409 done!\n",
      "g0047 done!\n",
      "g0028 done!\n",
      "g0325 done!\n",
      "g0417 done!\n",
      "g0154 done!\n",
      "g0122 done!\n",
      "g0254 done!\n",
      "g0021 done!\n",
      "g0314 done!\n",
      "g0248 done!\n",
      "g0187 done!\n",
      "g0227 done!\n",
      "g0164 done!\n",
      "g0115 done!\n",
      "g0223 done!\n"
     ]
    }
   ],
   "source": [
    "effe_lst=[];\n",
    "rat_flag_lst=[]\n",
    "for i in group_lst:  \n",
    "    nmi=pth_dh_xnhl+'/'+i+'_dh_coreg_pene.tif'\n",
    "    nm_gla=pth_gla+'\\\\'+i+'_gla.shp'\n",
    "    nm_dem_fill=pth_srtm_fill+'\\\\'+i+'_srtm.tif'\n",
    "    ds_dem=xr.open_dataset(nm_dem_fill)\n",
    "    ds_dh=xr.open_dataset(nmi,engine='rasterio')\n",
    "    df_gla=gu.Vector(nm_gla)\n",
    "    df_error=pd.read_csv(pth+'/xnhl_test/error_glacier/error_%s.csv'%i)\n",
    "    for gg_code in df_gla['Glc_Code']:\n",
    "        gg_shp = gu.Vector(df_gla.ds[df_gla.ds[\"Glc_Code\"] == gg_code])\n",
    "        gg_h=ds_dem.rio.clip(gg_shp.geometry)['band_data']\n",
    "        gg_dh0=ds_dh.rio.clip(gg_shp.geometry)['band_data']\n",
    "        gg_dh=gg_dh0.interp(x=gg_h.x.values,y=gg_h.y.values)\n",
    "        gg_dh_rm,flag,ratt=outlier_rm_ip(gg_dh,gg_h)\n",
    "        ####################################\n",
    "        gla_ma=gg_h.median().values  #物质平衡线高度，RGI中有些是错误的，因此，这里直接计算\n",
    "        mask=gg_h.values*0+1  #冰川区mask  \n",
    "        mask_acc=(gg_h.where(gg_h>=gla_ma)).values*0+1###积累区mask\n",
    "        mask_abla=(gg_h.where(gg_h<gla_ma)).values*0+1 ##消融区mask\n",
    "        acc_rt=(gg_dh_rm*mask_acc).count().values/(gg_h*mask_acc).count().values\n",
    "        abla_rt=(gg_dh_rm*mask_abla).count().values/(gg_h*mask_abla).count().values\n",
    "        if (acc_rt>=0.4)&(abla_rt>=0.5):\n",
    "            #有效冰川  可以进行插值处理\n",
    "            rtt_flag=1\n",
    "            ############################\n",
    "            (gg_h.fillna(-9999)).rio.to_raster('temp.tif')\n",
    "            ref_h=xdem.DEM('temp.tif')\n",
    "            gg_shp_index=gg_shp.rasterize(ref_h)\n",
    "            signal = xdem.volume.get_regional_hypsometric_signal(ddem=gg_dh_rm.values,\n",
    "                                                                 ref_dem=ref_h.data,\n",
    "                                                                 glacier_index_map=gg_shp_index,)\n",
    "            ############\n",
    "            ddem_filled = xdem.volume.norm_regional_hypsometric_interpolation(\n",
    "                    voided_ddem=gg_dh_rm, ref_dem=ref_h, glacier_index_map=gg_shp_index, \n",
    "                    regional_signal=signal.dropna())  \n",
    "            ########################\n",
    "            gg_dh_rm_fill=gg_dh_rm*1\n",
    "            gg_dh_rm_fill.values=[ddem_filled]\n",
    "\n",
    "            df_error_gla=df_error[df_error['gla_id']==gg_code]\n",
    "            outnm_nc=pth+'/xnhl_test/gla_level_dh_final/nc/'+gg_code+'.nc'\n",
    "            #outnm_tif=pth+'/xnhl_test/gla_level_dh_final/nc_tif/'+gg_code+'.tif'\n",
    "            try:\n",
    "                nc_out_stp1(gg_dh_rm,gg_dh_rm_fill,gg_h,i,df_error_gla['effective samples'].values,df_error_gla['dh_error'].values,outnm_nc)\n",
    "            except:\n",
    "                nc_out_stp1(gg_dh_rm,gg_dh_rm_fill,gg_h,i,df_error_gla['effective samples'].values[0],df_error_gla['dh_error'].values[0],outnm_nc)\n",
    "            ###########################################\n",
    "            effe_lst.append([i,gg_code,gg_dh_rm.mean().values,gg_dh_rm.median().values,\n",
    "                             (gg_dh_rm*mask_acc).mean().values,(gg_dh_rm*mask_abla).mean().values,\n",
    "                             gg_dh_rm_fill.mean().values,gg_dh_rm_fill.median().values,\n",
    "                             (gg_dh_rm_fill*mask_acc).mean().values,(gg_dh_rm_fill*mask_abla).mean().values,\n",
    "                            ])\n",
    "            \n",
    "        else:\n",
    "            rtt_flag=0\n",
    "            ##无效冰川，不生产dh\n",
    "        rat_flag_lst.append([i,gg_code,ratt,acc_rt,abla_rt,rtt_flag])  \n",
    "        \n",
    "    print(i,'done!')\n",
    "    #break"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 237,
   "id": "f51b3507-a908-4c5e-9476-1e65bc1eca5a",
   "metadata": {},
   "outputs": [],
   "source": [
    "dff_effe=pd.DataFrame(effe_lst,\n",
    "                      columns=['group_id','gla_id','mean','median','mean_acc','mean_abla',\n",
    "                               'mean_fill','median_fill','mean_fill_acc','mean_fill_abla'])\n",
    "effe_nm=pth+'/xnhl_test/dh-effe-sta_xnhl.csv'\n",
    "dff_effe.to_csv(effe_nm,)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 238,
   "id": "c742fea6-d8df-4f07-855e-7d61756d5246",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>group_id</th>\n",
       "      <th>gla_id</th>\n",
       "      <th>mean</th>\n",
       "      <th>median</th>\n",
       "      <th>mean_acc</th>\n",
       "      <th>mean_abla</th>\n",
       "      <th>mean_fill</th>\n",
       "      <th>median_fill</th>\n",
       "      <th>mean_fill_acc</th>\n",
       "      <th>mean_fill_abla</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101096E33298N</td>\n",
       "      <td>-2.2477363664235055</td>\n",
       "      <td>-4.501034729049512</td>\n",
       "      <td>-6.477638160943147</td>\n",
       "      <td>2.094463705826747</td>\n",
       "      <td>-2.532606404069561</td>\n",
       "      <td>-4.124494294405553</td>\n",
       "      <td>-6.559231598021659</td>\n",
       "      <td>1.5191851973447386</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101108E33292N</td>\n",
       "      <td>1.7055120859696382</td>\n",
       "      <td>1.4564942258694655</td>\n",
       "      <td>2.391291888114324</td>\n",
       "      <td>1.3545383365419552</td>\n",
       "      <td>1.7664486632688579</td>\n",
       "      <td>1.8403107749920027</td>\n",
       "      <td>2.213537529938661</td>\n",
       "      <td>1.3144241329375668</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101093E33294N</td>\n",
       "      <td>-6.184100287433032</td>\n",
       "      <td>-3.108037407906206</td>\n",
       "      <td>-8.328385448673215</td>\n",
       "      <td>-4.013504878815669</td>\n",
       "      <td>-5.655661475976183</td>\n",
       "      <td>-3.263412033998384</td>\n",
       "      <td>-7.636487289328092</td>\n",
       "      <td>-3.656494682871014</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101093E33290N</td>\n",
       "      <td>-12.320571940108392</td>\n",
       "      <td>-8.514097458653362</td>\n",
       "      <td>-13.887114096194361</td>\n",
       "      <td>-10.910683999631022</td>\n",
       "      <td>-11.220769667891624</td>\n",
       "      <td>-9.17825530459772</td>\n",
       "      <td>-11.557833082840348</td>\n",
       "      <td>-10.849738621979077</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101093E33286N</td>\n",
       "      <td>-5.804584990834597</td>\n",
       "      <td>-6.774280173105558</td>\n",
       "      <td>-4.26381972249485</td>\n",
       "      <td>-7.105395012465693</td>\n",
       "      <td>-5.1462419181403005</td>\n",
       "      <td>-5.558762480361313</td>\n",
       "      <td>-3.311638526646453</td>\n",
       "      <td>-6.99367470398026</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14423</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087667E27839N</td>\n",
       "      <td>-14.305744244040415</td>\n",
       "      <td>-13.877755675423119</td>\n",
       "      <td>-14.621832453321016</td>\n",
       "      <td>-13.964569668943893</td>\n",
       "      <td>-14.547323091119024</td>\n",
       "      <td>-13.791156150933388</td>\n",
       "      <td>-14.600766258843318</td>\n",
       "      <td>-14.493696269897402</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14424</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087727E27847N</td>\n",
       "      <td>-7.46914588281632</td>\n",
       "      <td>-9.174980926666477</td>\n",
       "      <td>-1.908601579676263</td>\n",
       "      <td>-11.397401374389457</td>\n",
       "      <td>-6.562540148286721</td>\n",
       "      <td>-8.538949697103124</td>\n",
       "      <td>-1.7660597981723278</td>\n",
       "      <td>-11.385741837398966</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14425</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087725E27842N</td>\n",
       "      <td>-7.190112791233111</td>\n",
       "      <td>-10.623819344207073</td>\n",
       "      <td>-1.009441812480414</td>\n",
       "      <td>-10.311663790603161</td>\n",
       "      <td>-5.994087978294833</td>\n",
       "      <td>-7.356001129410733</td>\n",
       "      <td>-1.4136050454047535</td>\n",
       "      <td>-10.75490488492074</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14426</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087663E27836N</td>\n",
       "      <td>-5.2157884990193475</td>\n",
       "      <td>-4.1854104648007375</td>\n",
       "      <td>-5.743342886753424</td>\n",
       "      <td>-4.643015163765209</td>\n",
       "      <td>-7.52223914079417</td>\n",
       "      <td>-9.43930445332672</td>\n",
       "      <td>-7.5348090814424</td>\n",
       "      <td>-7.50966920014594</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14427</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087741E27837N</td>\n",
       "      <td>1.5278226475482433</td>\n",
       "      <td>4.047111306130491</td>\n",
       "      <td>6.873917654954094</td>\n",
       "      <td>-3.690857980336514</td>\n",
       "      <td>1.7077992754278637</td>\n",
       "      <td>4.2879104101051855</td>\n",
       "      <td>7.351684564190758</td>\n",
       "      <td>-3.999619431219792</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>14428 rows × 10 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "      group_id          gla_id                 mean               median  \\\n",
       "0        g0018  G101096E33298N  -2.2477363664235055   -4.501034729049512   \n",
       "1        g0018  G101108E33292N   1.7055120859696382   1.4564942258694655   \n",
       "2        g0018  G101093E33294N   -6.184100287433032   -3.108037407906206   \n",
       "3        g0018  G101093E33290N  -12.320571940108392   -8.514097458653362   \n",
       "4        g0018  G101093E33286N   -5.804584990834597   -6.774280173105558   \n",
       "...        ...             ...                  ...                  ...   \n",
       "14423    g0341  G087667E27839N  -14.305744244040415  -13.877755675423119   \n",
       "14424    g0341  G087727E27847N    -7.46914588281632   -9.174980926666477   \n",
       "14425    g0341  G087725E27842N   -7.190112791233111  -10.623819344207073   \n",
       "14426    g0341  G087663E27836N  -5.2157884990193475  -4.1854104648007375   \n",
       "14427    g0341  G087741E27837N   1.5278226475482433    4.047111306130491   \n",
       "\n",
       "                  mean_acc            mean_abla            mean_fill  \\\n",
       "0       -6.477638160943147    2.094463705826747   -2.532606404069561   \n",
       "1        2.391291888114324   1.3545383365419552   1.7664486632688579   \n",
       "2       -8.328385448673215   -4.013504878815669   -5.655661475976183   \n",
       "3      -13.887114096194361  -10.910683999631022  -11.220769667891624   \n",
       "4        -4.26381972249485   -7.105395012465693  -5.1462419181403005   \n",
       "...                    ...                  ...                  ...   \n",
       "14423  -14.621832453321016  -13.964569668943893  -14.547323091119024   \n",
       "14424   -1.908601579676263  -11.397401374389457   -6.562540148286721   \n",
       "14425   -1.009441812480414  -10.311663790603161   -5.994087978294833   \n",
       "14426   -5.743342886753424   -4.643015163765209    -7.52223914079417   \n",
       "14427    6.873917654954094   -3.690857980336514   1.7077992754278637   \n",
       "\n",
       "               median_fill        mean_fill_acc       mean_fill_abla  \n",
       "0       -4.124494294405553   -6.559231598021659   1.5191851973447386  \n",
       "1       1.8403107749920027    2.213537529938661   1.3144241329375668  \n",
       "2       -3.263412033998384   -7.636487289328092   -3.656494682871014  \n",
       "3        -9.17825530459772  -11.557833082840348  -10.849738621979077  \n",
       "4       -5.558762480361313   -3.311638526646453    -6.99367470398026  \n",
       "...                    ...                  ...                  ...  \n",
       "14423  -13.791156150933388  -14.600766258843318  -14.493696269897402  \n",
       "14424   -8.538949697103124  -1.7660597981723278  -11.385741837398966  \n",
       "14425   -7.356001129410733  -1.4136050454047535   -10.75490488492074  \n",
       "14426    -9.43930445332672     -7.5348090814424    -7.50966920014594  \n",
       "14427   4.2879104101051855    7.351684564190758   -3.999619431219792  \n",
       "\n",
       "[14428 rows x 10 columns]"
      ]
     },
     "execution_count": 238,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dff_effe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 240,
   "id": "9820f7fc-47a8-4abb-8b1b-ae55a29c6525",
   "metadata": {},
   "outputs": [],
   "source": [
    "dff_flag=pd.DataFrame(rat_flag_lst,\n",
    "                      columns=['group_id','gla_id','glacier_ratio','acc_ratio','abla_ratio','flag'])\n",
    "flag_nm=pth+'/xnhl_test/dh-glacier-flag.csv'\n",
    "dff_flag.to_csv(flag_nm,)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 242,
   "id": "d6dc6389-4340-4167-8396-12edccb19188",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>group_id</th>\n",
       "      <th>gla_id</th>\n",
       "      <th>glacier_ratio</th>\n",
       "      <th>acc_ratio</th>\n",
       "      <th>abla_ratio</th>\n",
       "      <th>flag</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101096E33298N</td>\n",
       "      <td>0.713396</td>\n",
       "      <td>0.720497</td>\n",
       "      <td>0.706250</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101108E33292N</td>\n",
       "      <td>0.704026</td>\n",
       "      <td>0.474067</td>\n",
       "      <td>0.936523</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101093E33294N</td>\n",
       "      <td>0.755760</td>\n",
       "      <td>0.756881</td>\n",
       "      <td>0.754630</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101093E33290N</td>\n",
       "      <td>0.701107</td>\n",
       "      <td>0.633803</td>\n",
       "      <td>0.775194</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>g0018</td>\n",
       "      <td>G101093E33286N</td>\n",
       "      <td>0.783972</td>\n",
       "      <td>0.715278</td>\n",
       "      <td>0.853147</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17979</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087667E27839N</td>\n",
       "      <td>0.898799</td>\n",
       "      <td>0.931507</td>\n",
       "      <td>0.865979</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17980</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087727E27847N</td>\n",
       "      <td>0.734722</td>\n",
       "      <td>0.606648</td>\n",
       "      <td>0.863510</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17982</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087725E27842N</td>\n",
       "      <td>0.575290</td>\n",
       "      <td>0.378788</td>\n",
       "      <td>0.779528</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17983</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087663E27836N</td>\n",
       "      <td>0.760417</td>\n",
       "      <td>0.791667</td>\n",
       "      <td>0.729167</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17987</th>\n",
       "      <td>g0341</td>\n",
       "      <td>G087741E27837N</td>\n",
       "      <td>0.618781</td>\n",
       "      <td>0.607916</td>\n",
       "      <td>0.629769</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>14428 rows × 6 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "      group_id          gla_id  glacier_ratio  acc_ratio  abla_ratio  flag\n",
       "0        g0018  G101096E33298N       0.713396   0.720497    0.706250     1\n",
       "2        g0018  G101108E33292N       0.704026   0.474067    0.936523     1\n",
       "3        g0018  G101093E33294N       0.755760   0.756881    0.754630     1\n",
       "4        g0018  G101093E33290N       0.701107   0.633803    0.775194     1\n",
       "6        g0018  G101093E33286N       0.783972   0.715278    0.853147     1\n",
       "...        ...             ...            ...        ...         ...   ...\n",
       "17979    g0341  G087667E27839N       0.898799   0.931507    0.865979     1\n",
       "17980    g0341  G087727E27847N       0.734722   0.606648    0.863510     1\n",
       "17982    g0341  G087725E27842N       0.575290   0.378788    0.779528     1\n",
       "17983    g0341  G087663E27836N       0.760417   0.791667    0.729167     1\n",
       "17987    g0341  G087741E27837N       0.618781   0.607916    0.629769     1\n",
       "\n",
       "[14428 rows x 6 columns]"
      ]
     },
     "execution_count": 242,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dff_flag[(dff_flag.acc_ratio>=0.2)&(dff_flag.abla_ratio>=0.3)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3092b9d3-b8d0-4f51-adb7-9044161fc11d",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "836b820a-3f18-4ef2-a8ff-b069b881c64a",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "16b14232-263f-44f5-aa76-1a63a2c4e39b",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f1d425de-b725-4a5c-9dc5-71fe913e0f7f",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "py39",
   "language": "python",
   "name": "py39"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.19"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}