1# Library of McStas runtime perl functions 2# 3# This file is part of the McStas neutron ray-trace simulation package 4# Copyright (C) 1997-2006, All rights reserved 5# Risoe National Laborartory, Roskilde, Denmark 6# Institut Laue Langevin, Grenoble, France 7# 8# This program is free software; you can redistribute it and/or modify 9# it under the terms of the GNU General Public License as published by 10# the Free Software Foundation; version 2 of the License. 11# 12# This program is distributed in the hope that it will be useful, 13# but WITHOUT ANY WARRANTY; without even the implied warranty of 14# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15# GNU General Public License for more details. 16# 17# You should have received a copy of the GNU General Public License 18# along with this program; if not, write to the Free Software 19# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 20# 21use Config; 22 23use Math::Amoeba qw(MinimiseND); 24require "mccode_config.perl"; 25 26# Overload with user's personal config 27if ($ENV{"HOME"} && -e $ENV{"HOME"}."/.mcstas/mccode_config.perl") { 28 require $ENV{"HOME"}."/.mcstas/mccode_config.perl"; 29} 30 31our $optim_bestvalue=0; 32our $optim_variables; 33our @optim_datablock = (); 34our @optim_first; 35our @optim_youts=(); 36 37# must run optimization steps as scans in oder to extract Detector lines 38# from output: (mcrun::do_scan) mcrun -N 1 39sub minimize_function { 40 my (@p) = @_; # array of instrument parameters at optimization step 41 my $j; 42 # function must return a value to be maximized 43 $y = 0; 44 45 my @scanned = (); 46 my @minval = (); 47 my @maxval = (); 48 49 my $out = "$optim_iterations "; 50 my @youts=(); 51 52 # create a copy of the global $scan_info within limits 53 for($j = 0; $j < @{$scan_info->{VARS}}; $j++) { 54 if ($p[$j] > $scan_info->{MAX}[$j]) { $p[$j] = $scan_info->{MAX}[$j]; } 55 elsif ($p[$j] < $scan_info->{MIN}[$j]) { $p[$j] = $scan_info->{MIN}[$j]; } 56 $minval[$j] = $p[$j]; 57 $maxval[$j] = $p[$j]; 58 $scanned[$j] = $scan_info->{VARS}[$j]; 59 $out .= "$p[$j] "; # parameter values for this iteration step 60 } 61 62 # assemble $scan_info from optim_info and optimization step ($numpoints=1) 63 my $optim_info = { VARS => \@scanned, MIN => \@minval, MAX => \@maxval }; 64 65 # execute single iteration step (do_scan). returns [params_val I err ... I err] 66 ($datablock, $variables, @youts) = do_scan($optim_info); 67 my @vars = split(" ", $variables); 68 my @vals = split(" ", $datablock); 69 70 # search optim_names in variables 71 # get list of monitors to maximize and loop 72 my $found_monitor=0; 73 for ($j = @{$scan_info->{VARS}}; $j < @vars; $j += 2) { 74 if ($optim_iterations == 0) { $optim_first[$j] = 0; } 75 my $value = $vals[$j]; 76 if ($optim_first[$j] == 0 && $value != 0) { 77 $optim_first[$j] = abs($value); 78 } 79 if ($optim_first[$j] > 0) { 80 $value /= $optim_first[$j]; 81 } 82 if ($optim_flag > 1) { # all monitors 83 $y = $y + $value; # add all values for criteria 84 $found_monitor++; 85 } else { # selected monitors 86 my $i; 87 for($i = 0; $i < @optim_names; $i++) { 88 my $this_name = $optim_names[$i]; 89 # add each value of monitor to $y 90 if ($vars[$j] eq "$this_name" . "_I") { 91 $y = $y + $value; # add corresponding value to found name for criteria 92 $found_monitor++; 93 } 94 } 95 } 96 $out .= " $vals[$j] $vals[$j+1]"; 97 } # end for j 98 die "optimization: ERROR: selected component is not a monitor\n" unless $found_monitor; 99 $out .= " $y 0"; 100 push @optim_datablock, "$out\n"; 101 if ($optim_iterations == 0) { 102 $optim_variables = "Point " . $variables; 103 @optim_youts = @youts; 104 $optim_variables .= " criteria null"; 105 push @optim_youts, " (criteria,null)"; 106 } 107 if ($optim_iterations == 0 || $y >= $optim_bestvalue) { 108 @optim_best = @p; 109 $optim_bestvalue = $y; 110 } 111 $optim_iterations++; 112 113 # Handle writing of mcplot-compatible optim output file 114 $numpoints = $optim_iterations; 115 if (!$optfile) { 116 if ($MCSTAS::mcstas_config{'TEMP'} ne "no") { 117 require File::Temp; 118 ($OPT, $optfile) = File::Temp::tempfile("mcoptim"."_XXXX", SUFFIX => '.dat'); # throw file handle 119 } else { 120 $optfile = int(rand(9999)); 121 $optfile = "mcoptim"."_$optfile.dat"; 122 } 123 } else { 124 $OPT = new FileHandle; 125 } 126 open($OPT, ">$optfile"); 127 autoflush $OPT 1; 128 push @opt_out, "$optim_iterations $y 0\n"; 129 output_dat_header($OPT, "# ",$scan_info, '(criteria,null)', 'Point criteria null', $optfile); 130 print $OPT "@opt_out\n"; 131 close($OPT); 132 133 print "Optimization $optim_iterations criteria=$y\n"; 134 135 return -$y; 136} 137 138sub sighandler_optim { 139 my $signame = shift; 140 print "mcrun: Recieved signal $signame during optimization (iteration $optim_iterations).\n"; 141 if ($signame eq "INT" || $signame eq "TERM") { 142 save_optimization(0); 143 exit(0); 144 } else { 145 print "# Optimization history:\n"; 146 my $datablock = join(" ", @optim_datablock); # mcoptimlib-> @optim_datablock, @optim_youts, $optim_variable 147 print "# variables: $optim_variables\n"; 148 print "$datablock\n"; 149 } 150} 151 152# function to save optimization history 153sub save_optimization { 154 my ($final_call) = @_; 155 my ($j); 156 157 # constrain within limits 158 for($j = 0; $j < @{$scan_info->{VARS}}; $j++) { 159 if ($optim_best[$j] > $scan_info->{MAX}[$j]) { $optim_best[$j] = $scan_info->{MAX}[$j]; } 160 elsif ($optim_best[$j] < $scan_info->{MIN}[$j]) { $optim_best[$j] = $scan_info->{MIN}[$j]; } 161 } 162 163 # display result 164 print "Iteration $optim_iterations parameters:\n"; 165 166 for($j = 0; $j < @{$scan_info->{VARS}}; $j++) { 167 $i = $scan_info->{VARS}[$j]; # Index of variable to be scanned 168 print "$params[$i]=$optim_best[$j] "; 169 } 170 print "\n"; 171 172 # reactivate $data_dir and relaunch minimize_function to fill directory 173 # mcstas.sim is last iteration result 174 if ($data_dir_saved) { 175 $data_dir = $data_dir_saved; 176 print "Generating optimized configuration in --dir=$data_dir\n"; 177 minimize_function(@optim_best); # this also sets scan-info by calling do_scan 178 } 179 # output optimization results as a scan 180 $numpoints = $optim_iterations; 181 my $prefix = ""; 182 if($data_dir) { $prefix = "$data_dir/"; } 183 my $datfile = "mcstas.dat"; 184 my $simfile = $datfile; 185 $simfile =~ s/\.dat$//; # Strip any trailing ".dat" extension ... 186 $simfile .= $format_ext; # ... and add ".sim|m|sci" extension. 187 my $datablock = join(" ", @optim_datablock); # mcoptimlib-> @optim_datablock, @optim_youts, $optim_variable 188 189 # Only initialize / use $DAT datafile if format is PGPLOT 190 if ($MCSTAS::mcstas_config{'PLOTTER'} =~ /PGPLOT|McStas/i) { 191 my $DAT = new FileHandle; 192 open($DAT, ">${prefix}$datfile"); 193 output_dat_header($DAT, "# ", $scan_info, \@optim_youts, $optim_variables, $datfile); 194 print $DAT "$datablock\n"; 195 close($DAT); 196 } else { 197 $datfile = $simfile; # Any reference should be to the simfile 198 } 199 200 output_sim_file("${prefix}$simfile", $scan_info, \@optim_youts, $optim_variables, $datfile, $datablock); 201 202 print "Optimization file: '${prefix}$datfile'\nOptimized parameters: $optim_variables\n"; 203 undef($data_dir); 204 $numpoints = 1; 205} 206 2071; 208