SEd = input('What is the standard error of the contrast? ');
meand = input('What is the contrast mean? ');
nu = input('what are the degrees of freedom of the contrast? ');
 
likelihoodmax = 0;
theta(1) = meand - 5*SEd;
inc = SEd/100;
 
 for B = 1:1000
    theta(B) = theta(1) + (B-1)*inc;
    likelihood(B) = (1 +(meand - theta(B))^2/(nu*SEd^2))^(-(nu+1)/2);
   
    if likelihood(B) > likelihoodmax
        likelihoodmax = likelihood(B);
    end
 end
 for B = 1:1000
     likelihood(B) = likelihood(B)/likelihoodmax;
 end
 
 outofrange = meand - 6*SEd;
 
 begin8 = outofrange;
 begin32 = outofrange;
 end8 = outofrange;
 end32 = outofrange;
 for B = 1:1000
     if begin8 == outofrange
         if likelihood(B) > 1/8
             begin8 = theta(B);
         end
     end
     if begin32 == outofrange
         if likelihood(B) > 1/32
             begin32 = theta(B);
         end
     end
     if and(begin8 ~= outofrange, end8 == outofrange)
         if likelihood(B) < 1/8
             end8 = theta(B);
         end
     end
     if and(begin32 ~= outofrange, end32 == outofrange)
         if likelihood(B) < 1/32
             end32 = theta(B);
         end
     end
 end
 
 
  theta1 = input('What is the population contrast mean assumed by the first hypothesis? ');
  theta2 = input('What is the population contrast mean assumed by the second hypothesis? ');
  B1 = int16((theta1 - theta(1))/inc  + 1);
  B2 = int16((theta2 - theta(1))/inc  + 1);
  likelihoodratio  = likelihood(B1)/likelihood(B2)
    
     begin32
     end32
     begin8
     end8