/usr/share/psychtoolbox-3/PsychRadiometric/PsychISO2007MPE/ISO2007MPEBasicTest.m is in psychtoolbox-3-common 3.0.11.20131230.dfsg1-1build1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 | % ISO2007MPEBasicTest
%
% ****************************************************************************
% IMPORTANT: Before using the ISO2007MPE routines, please see the notes on usage
% and responsibility in PsychISO2007MPE/Contents.m (type "help PsychISO2007MPE"
% at the Matlab prompt.
% ****************************************************************************
%
% Test code for our implementation of the ISO 2007 broadband MPE standard.
%
% We don't have any known test cases to check, particularly since we can
% only measure light in the visible. Here we verify that a bright sunlight
% measured in Philly is below the MPE, but not by all that much. This is
% broadly consistent with what we find when we compare the same light to the ANSI
% standard for laser light (with a few assumptions to apply that standard to
% broadband light.)
%
% 6/26/13 dhb Wrote it.
%% Clear and close
clear; close all;
%% Wavelength sampling.
%
% Easiest to cover the whole range covered by the standard,
% that way we don't have to think.
S = [200 1 1301];
%% Load CIE functions.
load T_xyz1931
T_xyz = SplineCmf(S_xyz1931,683*T_xyz1931,S);
%% Load in a test spectrum
%
% This is a bright sunlight measured through a window and off
% of a white piece of paper in Philly.
% We only have measurements between 380 and 780 nm. In this
% example, we zero extend when we spline to the whole range.
load spd_phillybright
spd_phillybright = SplineSpd(S_phillybright,spd_phillybright,S,0);
photopicLuminancePhillyBrightCdM2 = T_xyz(2,:)*spd_phillybright;
PLOT_SPECTRUM = 0;
if (PLOT_SPECTRUM)
figure;
plot(SToWls(S),spd_phillybright,'r','LineWidth',2);
xlabel('Wavelength (nm)');
ylabel('Radiance (Watts/[sr-m2-wlinterval]');
end
%% Specify stimulus parameters
stimulusDiameterDegrees = 10;
stimulusAreaDegrees2 = pi*((stimulusDiameterDegrees/2)^2);
stimulusDurationSecs = 60*60;
%% Eye length, for conversion from radiance to retinal irradiance
%
% This will default to 17 if you don't pass it to the compute routine.
eyeLengthMm = 17;
%% Plot of weighting functions
[wls,weightingR,weightingA,weightingS,wls_R,rawWeightingR,wls_A,rawWeightingA,wls_S,rawWeightingS] = ISO2007MPEGetWeighings(S);
figure; clf; set(gcf,'Position',[400 500 1400 550]);
subplot(1,3,1); hold on
plot(wls_R,rawWeightingR,'r','LineWidth',3);
plot(wls,weightingR,'b:','LineWidth',2);
xlabel('Wavelength (nm)')
ylabel('R_lambda')
title('Thermal Hazard Function');
xlim([200 1500]);
subplot(1,3,2); hold on
plot(wls_A,rawWeightingA,'r','LineWidth',3);
plot(wls,weightingA,'b:','LineWidth',2);
xlabel('Wavelength (nm)')
ylabel('A_lambda')
title('Aphakic Photochemical Hazard Function');
xlim([200 1500]);
subplot(1,3,3); hold on
plot(wls_S,rawWeightingS,'r','LineWidth',3);
plot(wls,weightingS,'b:','LineWidth',2);
xlabel('Wavelength (nm)')
ylabel('S_lambda')
title('UV Radiation Hazard Function');
xlim([200 1500]);
%% Do the computations
fprintf('\nRunning tests for a sunlight measured in Philadelphia\n\n');
[IsOverLimit,ISO2007MPEStruct] = ISO2007MPECheckType1ContinuousRadiance(S,spd_phillybright,stimulusDurationSecs,stimulusAreaDegrees2,eyeLengthMm);
ISO2007MPEPrintAnalysis(IsOverLimit,ISO2007MPEStruct)
%% Anterior segment limit for convergent beams
% This just makes sure the routine properly throws an error.
try
fprintf('\n\nTesting error catch for anterior segment limit for convergent beams\n');
[ISO2007MPEStruct.val8_UWattsPerCm2,limit8_UWattsPerCm2] = ISO2007MPEComputeType1ContinuousAntConvrgUnweightedValue(S,NaN,stimulusDurationSecs);
catch
end
|