function result = th2(q)

%*****************************************************************************80
%
%% th2()
%
    q2 = q*q;
    q2n = 1.0;
    p = 1.0;
    while ( true )
        q2n = q2n*q2;
        if ( ( 1.0 + q2n ) == 1.0 )
            break
        end
        p = p*((1  - q2n)*(1 + q2n))*(1 + q2n);
    end
    result = 2*sqrt(sqrt(q))*p;

  return
end
