lib/model/shrinkage_module_subs.pm

   1 # {{{ include
   2
   3 start include statements
   4
   5 use Data::Dumper;
   6
   7 end include statements
   8
   9 # }}} include
  10
  11 # {{{ enable
  12
  13 start enable
  14
  15 $self -> {'enabled'} = 1;
  16
  17 end enable
  18
  19 # }}} enable
  20
  21 # {{{ disable
  22
  23 start disable
  24
  25 $self -> {'enabled'} = 0;
  26
  27 end disable
  28
  29 # }}} disable
  30
  31 # {{{ format_shrinkage_tables
  32
  33 start format_shrinkage_tables
  34
  35 my $omegas = $self -> {'problem'} -> nomegas;
  36
  37 my $eta_str = 'ID';
  38 my $eps_str = 'ID IWRES EVID';
  39
  40 for( my $j = 1; $j <= $omegas; $j++ ) {
  41   $eta_str = $eta_str.' ETA'.$j;
  42 }
  43 $eta_str = $eta_str.' FILE='.$self -> eta_tablename.
  44     ' NOAPPEND ONEHEADER NOPRINT FIRSTONLY';
  45 $eps_str = $eps_str.' FILE='.$self -> wres_tablename.
  46     ' NOAPPEND ONEHEADER NOPRINT'."\n";
  47
  48 my $eta_table = model::problem::table -> new ( record_arr => [$eta_str] );
  49 my $eps_table = model::problem::table -> new ( record_arr => [$eps_str] );
  50
  51 my $wrapped = ( defined $self -> {'model'} -> datas -> [$self -> problem_number-1] -> wrap_column and
  52                 defined $self -> {'model'} -> datas -> [$self -> problem_number-1] -> cont_column );
  53
  54 if( $wrapped ) {
  55   # No good contifying the eta table. We use FIRSTONLY which will make
  56   # CONT=1 on all rows and PsN does not like that.
  57   $eps_table -> contify;
  58 }
  59
  60 @formatted = ( @{$eta_table -> _format_record},
  61                @{$eps_table -> _format_record} );
  62
  63 end format_shrinkage_tables
  64
  65 # }}} format_shrinkage_tables
  66
  67 # {{{ eta_tablename
  68
  69 start eta_tablename
  70
  71 my $probnum = $self -> problem_number;
  72
  73 $filename = 'prob'.'_'.$probnum.'.psn_etas';
  74
  75 end eta_tablename
  76
  77 # }}} eta_tablename
  78
  79 # {{{ wres_tablename
  80
  81 start wres_tablename
  82
  83 my $probnum = $self -> problem_number;
  84
  85 $filename = 'prob'.'_'.$probnum.'.psn_wres';
  86
  87 end wres_tablename
  88
  89 # }}} wres_tablename
  90
  91 # {{{ eta_table_exists
  92
  93 start eta_table_exists
  94
  95 my $directory = $self -> {'model'} -> directory;
  96 my $filename = $self -> eta_tablename;
  97
  98 $exists = -e $directory.$filename ? 1 : 0;
  99
 100 end eta_table_exists
 101
 102 # }}} eta_table_exists
 103
 104 # {{{ wres_table_exists
 105
 106 start wres_table_exists
 107
 108 my $directory = $self -> {'model'} -> directory;
 109 my $filename = $self -> wres_tablename;
 110
 111 $exists = -e $directory.$filename ? 1 : 0;
 112
 113 end wres_table_exists
 114
 115 # }}} wres_table_exists
 116
 117 # {{{ problem_number
 118
 119 start problem_number
 120
 121 if( defined $self -> {'model'} -> problems ) {
 122   # This is the default
 123   my @modprobs = @{$self -> {'model'} -> problems};
 124   for( my $i = 0; $i <= $#modprobs; $i++ ) {
 125     $problem_number = ($i+1) if( $modprobs[$i] eq $self -> {'problem'} );
 126   }
 127 } else {
 128   # This happens when the problems are not yet set in the model
 129   $problem_number = $self -> {'temp_problem_number'};
 130 }
 131
 132 end problem_number
 133
 134 # }}} problem_number
 135
 136 # {{{ eta_shrinkage
 137
 138 start eta_shrinkage
 139
 140 # We do not handle subproblem eta shrinkage
 141
 142 if( $self -> {'enabled'} ) {
 143   my $omegas  = $self -> {'model'} -> outputs -> [0] -> omegas;
 144   my $probnum = $self -> problem_number - 1; # start index at 0
 145   my @omega_indexes;
 146   if( defined $self -> {'model'} -> outputs -> [0] and
 147       defined $self -> {'model'} -> outputs -> [0] -> problems -> [0] and
 148       defined $self -> {'model'} -> outputs -> [0] -> problems -> [0] -> omega_indexes ) {
 149     @omega_indexes = @{$self -> {'model'} ->
 150                            outputs -> [0] -> problems -> [0] -> omega_indexes};
 151   }
 152   my $defined_indexes = 0;
 153   foreach my $index ( @omega_indexes ) {
 154     $defined_indexes++ if defined $index;
 155   }
 156
 157   if( defined $omegas and defined $omegas -> [$probnum] ) {
 158     if( scalar @{$omegas -> [$probnum]} == 1 ) { # One subprob
 159       if( $defined_indexes and $self -> eta_table_exists ) {
 160         my $sh_table = data -> new( directory            => $self -> {'model'} -> directory,
 161                                     filename             => $self -> eta_tablename,
 162                                     ignore_missing_files => 1,
 163                                     target               => 'mem',
 164                                     table_file           => 1 );
 165         my $diag_omega_idx = 0;
 166         for( my $j = 0; $j < scalar @{$omegas -> [$probnum][0]}; $j++ ) {
 167           next if( $omega_indexes[$j][0] != $omega_indexes[$j][1] );
 168           if ( defined $omegas -> [$probnum][0][$j] and defined $sh_table and
 169                $omegas -> [$probnum][0][$j] != 0 ) {
 170             my $omega = sqrt($omegas -> [$probnum][0][$j]);
 171             my $eta_sd = $sh_table -> sd( column => ($diag_omega_idx+2) );
 172             my $shrinkage = ($omega - $eta_sd)/$omega;
 173             $eta_shrinkage[0][$diag_omega_idx] = $shrinkage;
 174           } else {
 175             $eta_shrinkage[0][$diag_omega_idx] = undef;
 176           }
 177           $diag_omega_idx++;
 178         }
 179       } else {
 180         $eta_shrinkage[0] = [];
 181       }
 182     } elsif( scalar @{$omegas -> [$probnum]} == 0 ) {
 183       debug -> die( message => $self -> {'model'} -> full_name. " Call to output->omegas is empty. PsN can not compute shrinkage." );
 184     } else {
 185       debug -> die( message => $self -> {'model'} -> full_name . " Call to output->omegas indicates that results ".
 186                     "exists in multiple subproblems.PsN can not yet compute shrinkage".
 187                     " on the subproblem level" );
 188     }
 189   }
 190 }
 191
 192 end eta_shrinkage
 193
 194 # }}} eta_shrinkage
 195
 196 # {{{ wres_shrinkage
 197
 198 start wres_shrinkage
 199
 200 # We do not handle subproblem wres shrinkage
 201
 202 if( $self -> {'enabled'} ) {
 203   my $ofv  = $self -> {'model'} -> outputs -> [0] -> ofv; # Use ofv to test success
 204
 205   my $probnum = $self -> problem_number - 1; # start index at 0
 206   if( defined $ofv ) {
 207     if( scalar @{$ofv -> [$probnum]} == 1 ) {
 208       my $sh_table;
 209       if( defined $ofv -> [$probnum][0] and $self -> wres_table_exists ) {
 210         $sh_table = data -> new( directory            => $self -> {'model'} -> directory,
 211                                  filename             => $self -> wres_tablename,
 212                                  ignore_missing_files => 1,
 213                                  target               => 'mem',
 214                                  table_file           => 1 );
 215       }
 216       if( defined $sh_table ) {
 217         my $wres_sd = $sh_table -> sd( column        => 2,
 218                                        subset_column => 3,
 219                                        subset_syntax => '==0',
 220                                        global_sd     => 1 );
 221         my $shrinkage = 1 - $wres_sd;
 222         $wres_shrinkage[0] = $shrinkage;
 223       } else {
 224         $wres_shrinkage[0] = undef;
 225       }
 226     } elsif ( @{$ofv -> [$probnum]} < 1 ) {
 227       debug -> warn( level => 1,
 228                      message => "There seems to be a problem with the results from ".
 229                      $self -> {'model'} -> filename().". Cannot compute shrinkage." );
 230     } else {
 231       debug -> die( message => "Call to output->omegas indicates that results ".
 232                     "exists in multiple subproblems.PsN can not yet compute shrinkage".
 233                     " on the subproblem level" );
 234     }
 235   }
 236 }
 237
 238 end wres_shrinkage
 239
 240 # }}} wres_shrinkage
 241
 242
 243
 244
 245
 246
 247