Elaboradar  0.1

◆ analyse_VPR()

int elaboradar::CalcoloVPR::analyse_VPR ( float *  vpr_liq,
int *  snow,
float *  hliq 
)

funzione che analizza il profilo

analizza il profilo usando : la temperatura al suolo, la quota del massimo, e una funzione di interpolazione

Parametri
[out]vpr_liqvalore del profilo al livello liquido
[out]snowmatrice che indica se è presente neve o no secondo l'analisi fatta
[out]hliqquota del livello liquido
Restituisce
ier_ana valore che indica se tutto è andato a buon fine (0) o no (1)

Definizione alla linea 1180 del file cum_bac.cpp.

1182 {
1183  int ier=1,ier_ana=0,liv0;
1184  char date[]="000000000000";
1185  struct tm *tempo;
1186  time_t Time;
1187 
1188  // ------------inizializzazione delle variabili ----------
1189 
1190  //strcpy(date,"000000000000");
1191 
1192  int tipo_profilo=-1;
1193  float v600sottobb=NODATAVPR;
1194  float v1000=NODATAVPR;
1195  float v1500=NODATAVPR;
1196  float vliq=NODATAVPR;
1197  float vhliquid=NODATAVPR;
1198  float vprmax=NODATAVPR;
1199  //*togliere gli ultimi tre*/;
1200 
1201  //ier_max=trovo_hvprmax(&hvprmax);
1202 
1203 
1204  if (t_ground == NODATAVPR) //1 se non ho nè T nè il massimo esco altrimenti tipo_profilo=0
1205  {
1206  LOG_WARN("non ho T,...");
1207 
1208  if ( ! ier_max ) {
1209  LOG_ERROR(" non ho trovato hvprmax, nè T, esco");
1210  return 1;
1211  }
1212  tipo_profilo=0;
1213  }
1214  else
1215  {
1216 
1217  if (t_ground >= T_MAX_ML+0.65*(float)(livmin+TCK_VPR/2)/100.){ //1 se T > T_MAX_ML e non ho il massimo esco
1218  if ( ! ier_max ) {
1219  LOG_ERROR(" temperatura alta e non ho trovato hvprmax, esco");
1220  return 1;
1221  }
1222  tipo_profilo=0;
1223  }
1224 
1225 
1226  // if (t_ground >= T_MAX_SN+0.65*(float)(livmin+TCK_VPR/2)/100 && t_ground < T_MAX_ML+0.65*(float)(livmin+TCK_VPR/2)/100. )
1227  if (t_ground >= T_MAX_SN && t_ground < T_MAX_ML+0.65*(float)(livmin+TCK_VPR/2)/100. )
1228  {
1229 
1230  if ( ier_max ) {
1231  LOG_INFO(" temperatura da scioglimento e massimo in quota");
1232  tipo_profilo=2;
1233  }
1234  else{
1235  LOG_ERROR(" temperatura da scioglimento ma superiore a temperatura max neve e non ho trovato hvprmax, esco");
1236  return 1;
1237  }
1238  // solo una scritta per descrivere cos'è accaduto
1239  liv0=livmin+HALF_BB;
1240  if (hvprmax > liv0) LOG_INFO(" il livello %i è sotto la Bright band, ma T bassa interpolo",livmin);
1241  else LOG_INFO(" il livello %i potrebbe essere dentro la Bright Band, interpolo",livmin);
1242 
1243  }
1244 
1245  //if (t_ground >= T_MIN_ML && t_ground < T_MAX_SN+0.65*(float)(livmin+TCK_VPR/2)/100.)
1246  if (t_ground < T_MAX_SN)
1247  {
1248  if ( ier_max ){
1249  LOG_INFO(" temperatura da neve o scioglimento e massimo in quota");
1250  tipo_profilo=2;
1251  }
1252  else {
1253  LOG_INFO(" temperatura da neve o scioglimento e massimo non trovato,neve , non interpolo");
1254  tipo_profilo=3;
1255  hvprmax=0;
1256  }
1257  }
1258 
1259  }
1260 
1261  // InterpolaVPR_NR iv;
1262  InterpolaVPR_GSL iv;
1263 
1264  switch
1265  (tipo_profilo)
1266  {
1267  case 0:
1268  case 1:
1269  case 2:
1270  ier=iv.interpola_VPR(vpr.val.data(), hvprmax, livmin);
1271  if (ier){
1272  LOG_INFO(" interpolazione fallita");
1273  switch (tipo_profilo)
1274  {
1275 
1276  /*Questo fallisce se la bright band non è al suolo (per evitare correzioni dannose in casi poco omogenei)*/
1277  case 0:
1278  case 1:
1279  ier_ana=1;
1280  *vpr_liq=NODATAVPR;
1281  *hliq=NODATAVPR;
1282  break;
1283  case 2:
1284  *vpr_liq=vpr.val[(hvprmax+1000)/TCK_VPR]*2.15;/*21 aprile 2008*/
1285  *hliq=0;
1286  break;
1287  }
1288  }
1289  else{
1290  LOG_INFO(" interpolazione eseguita con successo");
1291  //
1292  // stampa del profilo interpolato
1293  const char* vpr_arch = getenv("VPR_ARCH");
1294  if (!vpr_arch) throw runtime_error("VPR_ARCH is not defined");
1295  string fname(vpr_arch);
1296  fname += "_int";
1297  File file(logging_category);
1298  file.open(fname, "wt", "vpr interpolato");
1299  for (unsigned i = 0; i < NMAXLAYER; ++i)
1300  fprintf(file," %f \n", cum_bac.dbz.RtoDBZ(iv.vpr_int[i]));
1301 
1302  /*calcolo valore di riferimento di vpr_liq per l'acqua liquida nell'ipotesi che a[2]=quota_bright_band e a[2]-1.5*a[3]=quota acqua liquida*/
1303  if (tipo_profilo == 2 ) {
1304  *hliq=(iv.E-2.1*iv.G)*1000.;
1305  //lineargauss(a[2]-2.1*a[3], a, vpr_liq, dyda, ndata);
1306  if (*hliq<0)
1307  *hliq=0; /*con casi di bright band bassa.. cerco di correggere il più possibile*/
1308  *vpr_liq=vpr.val[(hvprmax+1000)/TCK_VPR]*2.15;
1309  }
1310  else {
1311  *hliq=(iv.E-2.1*iv.G)*1000.;
1312  //lineargauss(a[2]-2.1*a[3], a, vpr_liq, dyda, ndata);
1313  if ( *hliq > livmin) {
1314  *vpr_liq=vpr.val[(int)(*hliq/TCK_VPR)]; // ... SE HO IL VALORE VPR USO QUELLO.
1315  }
1316  else // altrimenti tengo il valore vpr neve + 6 dB* e metto tipo_profilo=2
1317  {
1318  if (*hliq<0) *hliq=0;
1319  tipo_profilo=2;
1320  //*vpr_liq=vpr.val[(hvprmax+1000)/TCK_VPR]*2.15;
1321  // iv.C = risultato dell'interpolazione (interpolated vpr)
1322  *vpr_liq=iv.C;
1323  }
1324  }
1325  }
1326  break;
1327  case 3:
1328  *snow=1;
1329  *vpr_liq=NODATAVPR;
1330  *hliq=NODATAVPR;
1331  break;
1332  }
1333  LOG_INFO("TIPO_PROFILO= %i vpr_liq %f hliq %f", tipo_profilo, *vpr_liq,*hliq );
1334 
1335 
1336  /* parte di stampa test vpr*/
1337 
1338  /* nome data */
1339  //definisco stringa data in modo predefinito
1340  Time = cum_bac.volume.load_info->acq_date;
1341  tempo = gmtime(&Time);
1342  sprintf(date,"%04d%02d%02d%02d%02d",tempo->tm_year+1900, tempo->tm_mon+1,
1343  tempo->tm_mday,tempo->tm_hour, tempo->tm_min);
1344  if (! ier ) {
1345  if(*hliq > livmin +200 )
1346  vhliquid=cum_bac.dbz.RtoDBZ(vpr.val[(int)(*hliq)/TCK_VPR]);
1347  vliq=cum_bac.dbz.RtoDBZ(*vpr_liq);
1348  }
1349  if (ier_max) {
1350  if ( hvprmax-600 >= livmin )
1351  v600sottobb=cum_bac.dbz.RtoDBZ(vpr.val[(hvprmax-600)/TCK_VPR]);
1352  if ((hvprmax+1000)/TCK_VPR < NMAXLAYER )
1353  v1000=cum_bac.dbz.RtoDBZ(vpr.val[(hvprmax+1000)/TCK_VPR]);
1354  if ((hvprmax+1500)/TCK_VPR < NMAXLAYER )
1355  v1500=cum_bac.dbz.RtoDBZ(vpr.val[(hvprmax+1500)/TCK_VPR]);
1356  vprmax=cum_bac.dbz.RtoDBZ(vpr.val[(hvprmax/TCK_VPR)]);
1357  }
1358 
1359  if (FILE* test_vpr=fopen(getenv("TEST_VPR"),"a+"))
1360  {
1361  fprintf(test_vpr,"%s %i %i -1 %f %f %f %f %f %f %f %f %f %f %f %f %f %f \n",date,hvprmax,tipo_profilo,iv.chisqfin,*hliq,vliq,vhliquid,v600sottobb,v1000+6,v1500+6,vprmax,iv.rmsefin,iv.B,iv.E,iv.G,iv.C,iv.F);
1362  fclose(test_vpr);
1363  }
1364 
1365  // fine parte di stampa test vpr
1366 
1367  //---SCRIVO ALTEZZA MASSIMO PER CLASSIFICAZIONE AL RUN SUCCESSIVO
1368 
1370  LOG_INFO("fatta scrittura hmax vpr = %d",hvprmax);
1371 
1372  return (ier_ana);
1373 }
void write_vpr_hmax(int hvprmax)
write in $VPR_HMAX the vpr peak's height.
Definition: assets.cpp:322
radarelab::Volume< double > & volume
Set to Z undetect value the Zpixels classified as non-meteo echoes.
Definition: cum_bac.h:106
Assets assets
others
Definition: cum_bac.h:88
radarelab::algo::DBZ dbz
????
Definition: cum_bac.h:110
Open a file taking its name from a given env variable.
Definition: utils.h:22
std::shared_ptr< LoadInfo > load_info
Polar volume information.
Definition: volume.h:272
int ier_max
flag d'errore su calcolo quota max
Definition: cum_bac.h:248
int livmin
quota livello minimo calcolato
Definition: cum_bac.h:243
radarelab::algo::VPR vpr
Informa se il pixel è convettivo.
Definition: cum_bac.h:236
int hvprmax
quota picco vpr
Definition: cum_bac.h:237
CUM_BAC & cum_bac
oggeto CUM_BAC di riferimento
Definition: cum_bac.h:230
float t_ground
2m temperature
Definition: cum_bac.h:233

Referenzia elaboradar::CUM_BAC::date, elaboradar::CUM_BAC::logging_category, e radarelab::File::open().