text:
# If in this problem I regress $x$ on $y$ instead than $y$ on $x$, do I need to use an error-in-variables model?

I was trying to write an answer for this question:

Selection of data range changes coefficients too much in lmer (inverse regression)

Basically the OP has lots of data of Amplification vs Voltage (see the image below) for different tested devices and would like to estimate, for each device, the Voltage given an Amplification = 150.

Apparently, in literature there's a physical model which has been used for fitting data from these devices:

$$\text{Amplification}=\frac{1}{1-\left(\frac{\text{Voltage}}{p_0}\right)^{p_1}}+p_2+\epsilon$$

However, the model has Amplification as the dependent variable and Voltage as the predictor. This means that if we fit a nonlinear mixed model to the data with nlmer (for example), then we will get confidence intervals on the Amplification given Voltage, not the other way around.

I thought of simply regressing Voltage on Amplification, instead, i.e., invert the above relationship to get

$$\text{Voltage}=f(\text{Amplification};\theta_1,\theta_2,\theta_3)+\epsilon = \theta_1\left(1-\frac{1}{\text{Amplification}+\theta_2}\right)^{\theta_3} +\epsilon$$

and fit a nonlinear mixed model with nlmer doing something like

model <- function(x, theta1, theta2, theta3) {
y <- theta1 * (1 - (1 / (x + theta2))) ^ theta3
}

starting_values <- c(theta1 = 1, theta2 = -180, theta3 = 0.8)
nlmer(APD_data, Voltage ~ model(Amplification, theta1, theta2, theta3) ~
theta1 + theta2 + theta3 + theta1|Serial_number + theta2|Serial_number +
theta3|Serial_number, start = starting_values, verbose = TRUE)

(Serial_number is the tag that identifies each APD device: think of it as a Subject variable).

However, the experiments were conducted (as far as I understood) by controlling Voltage and measuring Amplification. If I invert the two variables in the regression, does it mean I need to build an errors-in-variables model? In that case I give up - nonlinear mixed errors-in-variables models are way above my pay grade :)

• It depends on what kind of relation this 'voltage given an amplification = 150' is. Because of the errors you will get regression dilution (not only a problem for the direction of the relation), but this is not creating a bias if the task is to make a point prediction based on the estimated model (the measurement to which the point estimate applies has the same error and regression dilution). On the other hand, I wonder if the task is really to make those point estimates for the same situation, and the 'voltage at amp=150' may really be more like a sort of parameter estimate. – Sextus Empiricus Nov 15 '17 at 22:24
• @MartijnWeterings I'm not sure why regression dilution is not a problem for predicting from the estimated model... Amplification=150 is not necessarily a measurement. I'm also very interested by your comment about this prediction problem being really an estimation problem in disguise - do you mean one of the parameters in the original model may correspond to 'voltage at amp = 150'? For sure it's not the voltage where Amplification becomes infinite - that's larger. Anyway, if you're right, what would you suggest to do? Build a model which explicitly contains that voltage as a parameter? – DeltaIV Nov 16 '17 at 7:59
• Like done with growth curves that have halfway points or maximum points (or compare also with Michaelis-Menten kinetics), you could be able to parameterize those curves with some of it's points. For this type of curve it may not be easy to do it explicitly. Yet, it may be done implicitly, if the 150 amp voltage value is more used to characterize the device (for use elsewhere). Saying, 'this is a device with a 150 amp voltage of ...' is the same sort of like saying 'this is a device with a $p_1$ of ....'. – Sextus Empiricus Nov 16 '17 at 8:52
• The contrast is with using the 150 amp voltage value, $V_{150}$ , as a prediction for the amplification. The $E(amp|V_{150})$ will be 150 (if you put it in the same setting, that is same error in the measurements of Voltage and amplification, so it is not valid for giving the value to, say, a customer that buys the device and makes different errors ). – Sextus Empiricus Nov 16 '17 at 8:53
• The problem has remained around in my mind. I have to find time to cast it into a nicely formulated form, and also try some things out with the data of the specific problem. – Sextus Empiricus Nov 19 '17 at 13:06

### The direction of regression may be important to prevent attenuation bias.

Your question about the regression $x \sim y$ versus $y \sim x$ has many angles. A problem which you might encounter is regression attenuation or regression dilution.

• This does not depend on which variable the experiment was controlled in the experiment or in which direction the causal relation has.
• It does depend on the error made in the variables.

It happens if the 'independent variable in the regression' has a large error.

The underlying mathematics does not care about the direction of the causal relationship, or which variable was controlled for, it cares about the errors in the variables.

$$(y+\epsilon_y) = a + b(x + \epsilon_x) \qquad vs \qquad \frac{(y+\epsilon_y)-a}{b} = (x + \epsilon_x)$$

What we do when changing the "direction" of a regression, is not really changing the direction, but much more like ignoring either $\epsilon_x$ (situation left hand side) or $\epsilon_y$ (situation right hand side).

The variable that was controlled in the experiment, Voltage in this case, can also have a large error, even if it was "controlled". (you do not really set the voltage, you set some button or switch that controls the voltage and eventually you measure the voltage by reading it from a voltmeter or something). So this left hand side situation with the $\epsilon_x$ ignored $(y+\epsilon_y) = a + b(x)$ may be wrong and cause problems (that is attenuation).

### 'In this problem' you do not have to worry about direction of regression and attenuation.

In this problem there is not a large error, or at least not much noise. The curves are smooth with little jitter. If there is an error then it is a systematic error, but such errors are not linked to regression attenuation (which is about the random errors). Also such systematic errors have little to do with other issues in the mathematics. Except for possibly making some inverse regression ill-posed due to crossing some asymptote or creating negative values in roots logs etcetera.

These systematic errors are more like something that should be dealt with on the practical side (testing equipment, performing good calibration, etcetera).

### Much more important is to use the proper model.

1. In the referenced question I have shown how the polynomial model is not working well in the large range https://stats.stackexchange.com/a/315546/164061 .
2. The nonlinear model is not doing much better $$\text{Amplification}=\frac{1}{1-\left(\frac{\text{Voltage}}{p_0}\right)^{p_1}}+p_2+\epsilon$$

Or at least, I cant seem to make it converge. And I believe it is ill-posed.

3. All this work of trying this perfect fitting is a bit an overkill for the simple task of getting an estimate for the Voltage value at Amplification = 150. This is an interpolation problem, not a fitting problem! This can be done by getting the values of the nearest data points above and below the 150 and use the line between these points to make the estimate. There is no noise which makes this method work badly. If there would be noise then one could use a line trough several points. This interpolation working well, shows that the direction of the relationship is not really the issue here.

If one does wish to fit a reasonable curve then I believe it would be better to dig deeper into the mechanics of the system and use knowledge of the device to create a good curve fit, rather than some polynomial or simplistic literature curve, which are "just" experimental relationships which have little value for generalization and provide little information on the mechanics and inner workings of the devices (which may not be the prime goal of the experiments, but would be a bonus, and at least would increase the robustness of the fitting method).

It is hard to do this work by just gazing at the raw data without much information of the system, however, I was able to get a reasonable fit for a differential equation with two power terms.

$$\frac{\partial A}{\partial V} = a(A-1)^b + c(A-1)^d + \epsilon$$

This is a separable equation and can be solved both numerically (easy, e.g. by deSolve in R) or analytically (although this seems to involve the hypergeometric function $_2F_1$ which is not easy to fit directly).

### Example of fitting with a differentiated function

Upon request the code to fit with the differentiated function.

I want to stress two points:

• It is indeed an interesting method to fit, although it a whole different question and not necessary (overkill) to make the predictions at Amplification 150. If there is more noise then it may even be be worse than a local linear or polynomial fit in the region around Amplification =150. The fit is only empirically determined and we can not guarantee that it helps to cancel out noise without introducing more bias in return. Some more information about the system (ab initio approach) would certainly help.

• I started out with this method of differentiation, in the first place to investigate and explore potential relationships based on differential equations. To use this method as a final fitting procedure is not always successful and recommendable, because the differentiation will amplify the noise. So it is not a general method to solve the problem and each problem has it's own quirks with a different approach (this also makes the answer to the question about 'regress $x$ on $y$ instead $y$ on $x$' not general and difficult to formulate).

code:

#### demonstrating fit based on differentiated data or differential equation
####
####    note that this is no production code
####    so no checks all sorts and several hard coded limits

library(deSolve)
library(optimr)

# we make two plots side by side
layout(matrix(c(1,2), 1, 2, byrow = TRUE))

# getting data (for simplicity of the example we just take the 91200913 serial number)

A <- c(1.00252,1.00452,1.00537,1.0056,1.00683,1.0069,1.00847,1.00935,1.01157,1.01418,1.01914,1.0247,1.02919,1.03511,1.04545,1.07362,1.11549,1.17123,1.25019,1.36276,1.5104,1.69862,1.9518,2.26756,2.66278,3.14247,3.73163,4.46152,5.36262,6.49514,7.9227,9.73803,12.0663,15.0943,19.1004,20.0563,21.0672,22.142,23.2867,24.5037,25.8102,27.2024,28.6916,30.2968,32.0181,33.8775,35.8937,38.0569,40.4069,42.9713,45.7766,48.8312,52.2068,55.916,60.0356,64.6109,69.7152,75.4698,82.0003,89.4222,97.9493,107.807,119.441,133.2,149.796,170.113,195.89,229.058,273.481,335.96,431.682,593.091,918.112,1903.74)
V <- c(24.9681,29.9591,34.9494,44.9372,49.9329,54.9625,59.9639,64.9641,69.965,74.9663,79.969,84.9719,89.974,94.9752,99.9759,109.974,119.969,129.96,139.96,149.963,159.958,169.959,179.963,189.957,199.97,209.971,219.971,229.973,239.966,249.962,259.971,269.971,279.97,289.968,299.959,301.968,303.967,305.966,307.965,309.955,311.963,313.963,315.961,317.962,319.956,321.951,323.961,325.962,327.963,329.965,331.966,333.959,335.97,337.972,339.973,341.974,343.967,345.97,347.978,349.98,351.98,353.971,355.983,357.983,359.983,361.975,363.989,365.989,367.984,369.979,371.992,373.994,375.985,377.999)

# numerical differentiation (the noise in this example is not high, alternatively one could higher order methods e.g. a Savitky Golay filter)

dA <- A[-1]-A[-74]
dV <- V[-1]-V[-74]
mA <- 0.5*(A[-1]+A[-74])
mV <- 0.5*(V[-1]+V[-74])

# plotting the derivartive dA/dv with a function of (A-1)
#     substracting A by 1, resulting in A-1 was done
#     to get an asymptote to zero instead of one.
#     this give probably more interresting relations on a log-log scale

col <- hsv(0,0,0.5)
plot(mA-1, dA/dV,
ylim=c(0.000001,2000), xlim=c(0.001,4000),
pch=21,cex=0.5, col = col, bg = col,
log="xy",xlab="amplification",ylab="d amplification / d Volt")
title("dA/dV as function of A \n data points and fit")

# this above plot looks so much like two seperate straight lines
#       lets try to fit a double power law
#
#       in the nls fit we use a scaling by the amplification
#       to give more weight on the lower values
#       (this arbitrary scaling does introduce some subjectivity,
#       but it is just a practical way to avoid the alternatively
#       plotting of logarithms which require the use of the 'port'
#       algorithm and limits to prevent logs of negative values.
#       so We can do it better, but we ar just lazy for this demonstration)

fit <- nls(dA/dV ~ a*(mA-1)^b+c*(mA-1)^d,
start=c(a=0.026, b=0.9, c=0.0005, d=1.88),
weights = (dA/dV)^-1,
control = nls.control(maxiter = 2000 , minFactor = 10^-9))
# ploting

coefs <- coef(fit)
x <- 10^(c(-10:16)/3)
lines(x,coefs[1]*x^coefs[2]+coefs[3]*x^coefs[4],col=col)

text_string <- paste0("fitted line: dA/dV =",round(coefs[1],5),"A^",round(coefs[2],2)," + ",round(coefs[3],5),"A^",round(coefs[4],2))
text(10^-3,10^-6,text_string,pos = 4)

# going back to the case A vs V
#     - we will need to integrate the fit
#
#     - note that the fit of dA/dV vs A was done to explore the relationship
#       the differentiation is not a very robust operation (noisy) and usually
#       we would whish to obatain a symbolic integration of the differential equation
#       and use the result to relate A as a function of V rather than dA/dV as a function of A
#
#       but this involves a hypergemoetric function 2F1 which happens to be difficult to fit
#       also the function is not that noisy at all so the fit in the space (dA/dV, A) works well
#
#     - in the end we do add some computation with the optim library to optimize the fit in the
#       space A vs V, thus without performing differentiation
#       (or at least, under the hood the optim algorithm will do some differtiation of the Loss function
#       to determine the convergence, but we do not do the noisy differntiation dA/dV)

# plotting data points
plot(V, A,
xlim=c(0,400), ylim=c(1,4000), log="y",
pch=21,cex=0.7, col = col, bg = col,
xlab="Voltage",ylab="Amplification")
title("A as function of V \n data points and two fits")

legend(0,4000,
c("data points","fit dA/dV vs A","fit A vs V"),
col=c(8,4,2),
pt.bg = c(col,0,0),
lty=c(NA,2,2),
pch=c(21,NA,NA))

# differential equation for use with deSolve
fitje <- function(t, state, parameters) {
with(as.list(c(state, parameters)), {
dV <-  -1
dA <-  -(a*(A-1)^b+c*(A-1)^d)
list(c(dV, dA))
})
}

# integration using deSolve (we put this in a function for easier use with optim)
integraal <- function(coefs) {
parameters <- coefs
state      <- c(V = max(V), A = max(A))
times      <- seq(0, 1900, by = 0.1)

out <- ode(y = state, times = times, func = fitje, parms = parameters)
out
}

out <- integraal(coefs)

# plot deSolve result
lines(out[,"V"],out[,"A"],
col=4,lty=2)

# solving the fit in the space V,A by putting the deSolve function into optim
#         (this will take potentially take some time and is only recommended if we can not
#          find an analytical function to put into optim/nls/nlmer/etc)

f <- function(coefs) {
A_o <- A
model <-  integraal(coefs)
A_m <- approx(x = model[,"V"], y = model[,"A"], xout = V)$y sum((A_o-A_m)^2) } f(coefs) # fingers crossed that it won't diverge to NAs due to getting out of the range of the approx-function # # not that using the coefficient from the previous fit will help the convergence # performing this method sec without a good initial guess may end badly model2 <- optim(coefs, f) #plotting coefs2 <- model2$par
out2 <- integraal(coefs2)

lines(out2[,"V"],out2[,"A"],
col=2,lty=2)

Image generated by this code:

• Impressive answer (+1). I will study it in detail, but for starters, why do you think the nonlinear model is ill-posed? – DeltaIV Nov 25 '17 at 20:49
• it doesn't fit (at least i cant get it to fit) – Sextus Empiricus Nov 25 '17 at 21:32
• The data clearly indicate asymptotes at amplifiction=1 and Voltage ~ 380. This fixes the coefficient $p_0$ and $p_2$, leaving only the power term $p_1$ to determine the rest of the curve, which does not succeed. Then those $p_0$ and $p_2$, which should have the values of the asymptotes, are gonna be changed in order to make the function still somewhat fitting. The best you get is a reasonable fit, but with the underlying principle (some realistic function) being violated because the fit only makes sense as an empirical relation and not theoretic. More likely is no fit due to bad convergence. – Sextus Empiricus Nov 25 '17 at 23:25
• Voltage ~ 380 is the industrial current voltage... interesting: I wonder if it might be related to the experimental setup. Amplification ~ 1 instead is intuitive - when there's no voltage applied, no amplification results. One last clarification - how do you fit the differential equation with the noise term? Could you please include the code you used to make the fit? Since $a,b,c,d$ are unknown, I think you need to estimate them while numerically integrating the ODE at the same time... – DeltaIV Nov 26 '17 at 14:24
• I have added the code. I changed one thing in the fit which is the weights which were previously defined by a function of the domain A^-2.5 and are now defined by a function of the image (dA/dV)^-1. – Sextus Empiricus Nov 28 '17 at 13:43


text:
# Differential Equations

 Question 1
Consider the following differential equation
$(1+y)\frac{dy}{dx}=y$
The solution of the equation that satisfies the condition is $y(1)=1$ is
 A $2ye^y=e^x+e$ B $y^2e^y=e^x$ C $ye^y=e^x$ D $(1+y)e^y=2e^x$
GATE ME 2021 SET-2   Engineering Mathematics
Question 1 Explanation:
\begin{aligned} (1+y) \frac{d y}{d x} &=y \\ \Rightarrow\qquad \left(\frac{1}{y}+1\right) d y &=d x \\ \Rightarrow\qquad \log y+y &=x+c \\ \text { Using, } \qquad y(1) &=1 \\ \quad \log 1+1 &=1+c \quad \Rightarrow c=0 \\ \text { Hence, } \quad \log y+y &=x \\ \Rightarrow\qquad \log y+y \operatorname{loge} &=x \\ \log _{\mathrm{e}}\left(y \cdot e^{y}\right) &=x \\ \Rightarrow\qquad y e^{y} &=e^{x} \end{aligned}
 Question 2
If the Laplace transform of a function $f(t)$ is given by $\frac{s+3}{(s+1)(s+2)}$ , then $f(0)$ is
 A 0 B $\frac{1}{2}$ C 1 D $\frac{3}{2}$
GATE ME 2021 SET-2   Engineering Mathematics
Question 2 Explanation:
By using partial fraction concept.
\begin{aligned} f(t) &=L^{-1}\left[\frac{s+3}{(s+1)(s+2)}\right] \\ &=L^{-1}\left[\frac{2}{s+1}-\frac{1}{s+2}\right] \\ \Rightarrow \qquad f(t) &=2 e^{-t}-e^{-2 t} \\ \text { So, } \qquad f(c)&=2 e^{0}-e^{0}=2-1=1 \end{aligned}
 Question 3
The ordinary differential equation $\frac{dy}{dt}=-\pi y$ subject to an initial condition $y(0)=1$ is solved numerically using the following scheme:

$\frac{y(t_{n+1})-y(t_n)}{h}=-\pi y(t_n)$

where $h$ is the time step, $t_n=nh,$ and $n=0,1,2,...$. This numerical scheme is stable for all values of $h$ in the interval.
 A $0 \lt h \lt \frac{2}{\pi}$ B $0 \lt h \lt 1$ C $0 \lt h \lt \frac{\pi}{2}$ D for all $h \gt 0$
GATE ME 2021 SET-1   Engineering Mathematics
Question 3 Explanation:
\begin{aligned} \frac{y\left(t_{n+1}\right)-y\left(t_{n}\right)}{h} &=-\pi y\left(t_{n}\right) \\ y_{n+1} &=-\pi / y_{n}+y_{n}=(-\pi h+1) y_{n} \end{aligned}
It is recursion relation between $y_{n+1}$ and $y_{n}$
So solution will be stable if
\begin{aligned} |-\pi h+1| & \lt 1 \\ -1 \lt -\pi h+1 & \lt 1 \\ -2 \lt -\pi h & \lt 0 \\ 0 & \lt \pi h \lt 2 \\ 0 & \lt h \lt \frac{2}{\pi} \end{aligned}
Therefore option (A) is correct.
 Question 4
The Dirac-delta function $(\delta (t-t_0)) \text{ for }t,t_0 \in \mathbb{R},$ has the following property

$\int_{a}^{b}\varphi (t)\delta (t-t_0)dt=\left\{\begin{matrix} \varphi (t_0) & a \lt t_0 \lt b\\ 0 &\text{otherwise} \end{matrix}\right.$

The Laplace transform of the Dirac-delta function $\delta (t-a)$ for $a \gt 0$;
$\mathcal{L} (\delta (t-a))=F(s)$ is
 A 0 B $\infty$ C $e^{sa}$ D $e^{-sa}$
GATE ME 2021 SET-1   Engineering Mathematics
Question 4 Explanation:
\begin{aligned} \because \qquad \int_{0}^{-} f(t) \delta(t-a) d t&=f(a) \\ \therefore \qquad L\{\delta(t-a)\}&=\int_{0}^{-} e^{-s t} \delta(t-a) d t=e^{-a s} \end{aligned}
 Question 5
If $y(x)$ satisfies the differential equation

$(\sin x)\frac{dy}{dx}+y \cos x =1$

subject to the condition $y(\pi /2)=\pi /2$, then $y(\pi /6)$ is
 A 0 B $\frac{\pi}{6}$ C $\frac{\pi}{3}$ D $\frac{\pi}{2}$
GATE ME 2021 SET-1   Engineering Mathematics
Question 5 Explanation:
\begin{aligned} \frac{d y}{d x}+y \cot x&=\text{cosec} x\\ 1.F. \qquad&=e^{\int \cot x d x}=e^{\log \sin x}=\sin x\\ \Rightarrow \quad y(\sin x)&=\int \text{cosec} x \sin x d x+c\\ \Rightarrow \qquad y \sin x&=x+c\\ \Rightarrow \qquad \frac{\pi}{2} \sin \frac{\pi}{2} & =\frac{\pi}{2}+c \\ \Rightarrow \qquad \frac{\pi}{2} & =\frac{\pi}{2}+c \quad \Rightarrow c=0 \\ \Rightarrow \qquad y \sin x & =x\\ \Rightarrow \qquad y \sin \frac{\pi}{6}&=\frac{\pi}{6}\\ \Rightarrow \qquad y\left(\frac{1}{2}\right) &=\frac{\pi}{6} \\ \Rightarrow y &=\frac{\pi}{3} \end{aligned}
 Question 6
The solution of
$\frac{d^2y}{dt^2}-y=1,$
which additionally satisfies $y|_{t=0}=\left.\begin{matrix} \frac{dy}{dt} \end{matrix}\right|_{t=0}=0$ in the Laplace s-domain is
 A $\frac{1}{s(s+1)(s-1)}$ B $\frac{1}{s(s+1)}$ C $\frac{1}{s(s-1)}$ D $\frac{1}{s-1}$
GATE ME 2020 SET-2   Engineering Mathematics
Question 6 Explanation:
\begin{aligned} y^{\prime\prime}-y &=1 \\ y(0) &=1 \\ y^{\prime}(0) &=1 \\ L\left\{y^{\prime\prime}-y\right\} &=L\{1\} \\ s^{2} Y(s)-s y(0)-y^{\prime}(0)-y(s) &=\frac{1}{s} \\ y(s) &=\frac{1}{s\left(s^{2}-1\right)}\\ & =\frac{1}{s(s+1)(s-1)} \end{aligned}
 Question 7
The Laplace transform of a function f(t) is $L(f)=\frac{1}{s^2+\omega ^2}$. Then f(t) is
 A $f(t)=\frac{1}{\omega ^2}(1-\cos \omega t)$ B $f(t)=\frac{1}{\omega} \cos \omega t$ C $f(t)=\frac{1}{\omega} \sin \omega t$ D $f(t)=\frac{1}{\omega ^2}(1-\sin \omega t)$
GATE ME 2020 SET-1   Engineering Mathematics
Question 7 Explanation:
$L(t)=\frac{1}{s^{2}+\omega ^{2}}$
$f(t)=L^{-1}\left \{ \frac{1}{s^{2}+\omega ^{2}} \right \}=\frac{1}{\omega }\sin \omega t$
 Question 8
A differential equation is given as
$x^2\frac{d^2y}{dx^2}-2x\frac{dy}{dx}+2y=4$
The solution of the differential equation in terms of arbitrary constants $C_1 \; and \; C_2$ is
 A $y=C_1 x^2+C_2 x+2$ B $y=\frac{C_1}{x^2}+C_2x+2$ C $y=C_1 x^2+C_2 x+2$ D $y=\frac{C_1}{x^2}+C_2x+4$
GATE ME 2019 SET-2   Engineering Mathematics
Question 8 Explanation:
Given $\left(x^{2} D^{2}-2 x D+2\right) y=4\dots(1)$
where $D=\frac{d}{d x}$
Let $x=e^{z}$ or $\log x=z \& x D=\theta, x^{2} D^{2}=\theta\left(\theta^{2}-1\right)\dots(2)$
Where $\theta=\frac{d}{d z}$
Using (2),(1) becomes
$\begin{array}{l} {[\theta(\theta-1)-2 \theta+2] y=4} \\ \Rightarrow\left(\theta^{2}-3 \theta+2\right) y=4 \\ \Rightarrow f(\theta) y=Q(z), \\ \text { where } f(\theta)=\theta^{2}-3 \theta+2 \& Q(z)=4 \end{array}$
C.F
Auxiliary equation if f(m)=0
$\begin{array}{l} \Rightarrow \mathrm{m}^{2}-3 \mathrm{m}+2=4 \\ \Rightarrow \mathrm{m}=1,2 \\ \Rightarrow \mathrm{y}_{\mathrm{C}}=\mathrm{C}_{1} \mathrm{e}^{2 \mathrm{z}}+\mathrm{C}_{2} \mathrm{e}^{\mathrm{z}} \\ \quad \mathrm{y}_{\mathrm{C}}=\mathrm{C}_{1} \mathrm{x}^{2}+\mathrm{C}_{2} \mathrm{x} \end{array}$
P.I:
$\theta(z)=4=4 \cdot e^{0 z}=k e^{a z+b}$
Here, $f(\theta)=f(a)=f(0)=0-0+2=2$
$\therefore y_{p}=\frac{1}{\int(a)} \theta(z)=\frac{1}{2} 4=2$
Hence, the solution of (1) is $y=y_{C}+y_{p}$
i.e., $y=\left(C_{1} x^{2}+C_{2} x\right)+2$
 Question 9
The differential equation $\frac{dy}{dx}+4y=5$ is valid in the domain $0\leq x\leq 1$ with y(0)=2.25. The solution of the differential equation is
 A $y=e^{-4x}+5$ B $y=e^{-4x}+1.25$ C $y=e^{4x}+5$ D $y=e^{4x}+1.25$
GATE ME 2019 SET-2   Engineering Mathematics
Question 9 Explanation:
$\begin{array}{l} \text { Given } \frac{\mathrm{dy}}{\mathrm{dx}}+4 \mathrm{y}=5,0 \leq \mathrm{x} \leq 1 \; \; \ldots(1) \\ \because \frac{d y}{d x}+P(x, y)=Q(x) \\ \text { With } \mathrm{y}(0)=2.25\;\;\; \ldots(2) \\ \text { Here, I.F }=\mathrm{e}^{\int 4 \mathrm{dx}}=\mathrm{e}^{4 \mathrm{x}} \end{array}$
The general solution of (1) is given by
$y \cdot e^{4 x}=\int(5)\left(e^{4 x}\right) d x+c$
$\Rightarrow \mathrm{y} \cdot \mathrm{e}^{4 \mathrm{x}}=\frac{5}{4} \mathrm{e}^{4 \mathrm{x}}+\mathrm{c}\cdots(3)$
Using (2) and (3)
$(2.25)(1)=\left(\frac{5}{4}\right)(1)+c$
$c=1 \cdots(4)$
The solution of (1) from (3) and (4) is
$y \cdot e^{4 x}=\frac{5}{4} e^{4 x}+1$
or $y=\frac{5}{4}+e^{4 x}=e^{-4 x}+1.25$
 Question 10
A harmonic function is analytic if it satisfies the Laplace equation.
If $u(x,y)=2x^2-2y^2+4xy$ is a harmonic function, then its conjugate harmonic function $v(x,y)$ is
 A $4xy-2x^2+2y^2+$ constant B $4y^2-4xy+$ constsnt C $2x^2-2y^2+xy+$ constsnt D $-4xy+2y^2-2x^2+$ constant
GATE ME 2019 SET-1   Engineering Mathematics
Question 10 Explanation:
\begin{aligned} &\text { Given } \mathrm{u}(\mathrm{x}, \mathrm{y})=2 \mathrm{x}^{2}-2 \mathrm{y}^{2}+4 \mathrm{xy}\\ &\frac{\partial u}{\partial x}=4 x+4 y, \frac{\partial u}{\partial y}=-4 y+4 x\\ &\text { By C-R equations } \frac{\partial u}{\partial x}=\frac{\partial v}{\partial y} \& \frac{\partial u}{\partial y}=\frac{-\partial v}{\partial x}\\ &\therefore \frac{\partial v}{\partial x}=\frac{-\partial u}{\partial y}=4 y-4 x \; \; \ldots(1)\\ &\frac{\partial v}{\partial y}=\frac{\partial u}{\partial x}=4 x+4 y \; \; \ldots(2)\\ &\begin{array}{l} \text { From }(1),\\ V=\int(4 y-4 x) d x+k(y) \\ V=4 x y-2 x^{2}+K(y) \;\; \ldots(3) \end{array}\\ &\therefore \frac{\partial v}{\partial y}=4 x+\frac{d}{d y} k(y)\ldots(4)\\ &\text { From }(2) \&(3), \\ &\frac{\mathrm{d}}{\mathrm{dy}} \mathrm{k}(\mathrm{y})=4 \mathrm{y}\\ &\int \frac{\mathrm{d}}{\mathrm{dy}} \mathrm{k}(\mathrm{y}) \mathrm{d} \mathrm{y}=\int 4 \mathrm{ydy}\\ &\mathrm{K}(\mathrm{y})=2 \mathrm{y}^{2}+\mathrm{C} \end{aligned}
while C is arbitrary real constant.
$\therefore$ Conjugate harmonic of u is
$\mathrm{V}=4 \mathrm{xy}-2 \mathrm{x}^{2}+2 \mathrm{y}^{2}+\mathrm{C}$
There are 10 questions to complete.

format: plain


text:
## NEC – vertical monopole radiation resistance study

A recent online discussion contained an analysis of the radiation efficiency of a vertical monopole over real ground.

The poster dismissed the values calculated by 4NEC2 and proposed his own formula $$RadiationEfficiency=\frac{35.6}{\mathbb{R}Z_f}$$ where 35.6 is the radiation resistance Rr of a quarter wave monopole over a perfectly conducting earth (PCE).

The reasoning seems to depend on Rr being independent of the ground type, but that is quite flawed.

## NEC insight

Let us look at an example of a quarter wave monopole with 120 shallow buried radials, average ground, at 3.8MHz.

Above is the model geometry. Continue reading NEC – vertical monopole radiation resistance study

## 4NEC2 – summary statistics discussion

I am a frequent user of 4NEC2 despite its many defects. It is a great work of software, in need of improvement that I suspect will never happen. Its author explained some time ago that it was developed in VB6 and with Windows upgrades, he no longer has a working VB6 development platform.

So, despite its defects, it is a very useful tool.

I will sometimes substitute _ for – in some 4NEC2 quantity labels for clarity in mathematical expressions. Continue reading 4NEC2 – summary statistics discussion

## FT240-43 matching transformer for an EFHW – NEC model at 3.6MHz

The article End Fed Half Wave matching transformer – 80-20m laid out a design for a EFHW transformer based on the readily available FT240-43.

This article builds an NEC model for an EFHW antenna at 3.6MHz incorporating a realistic model of the above transformer.

NEC provides for a NT card characterising a two port network using Y parameters.

## Y parameter model for the transformer

The Y parameter model is based on measured input impedance with port 2 open circuit, and short circuit, and the observed turns ratio.

Impedance was measured with the transformer at 3.6MHz using an AA-600.

Above, the calculated Y parameter model including a prototype NT card. This model captures the various loss components of the transformer, mainly magnetising loss, at 3.6MHz. Continue reading FT240-43 matching transformer for an EFHW – NEC model at 3.6MHz

## Small efficient matching transformer for an EFHW – NEC model at 3.6MHz

The article Small efficient matching transformer for an EFHW laid out a design for a small EFHW transformer.

This article builds an NEC model for an EFHW antenna at 3.6MHz incorporating a realistic model of the above transformer.

NEC provides for a NT card characterising a two port network using Y parameters.

## Y parameter model for the transformer

The Y parameter model is based on measured input impedance with port 2 open circuit, and short circuit, and the observed turns ratio.

Impedance was measured with the uncompensated transformer at 3.6MHz using an AA-600, the compensation in the reference article has little effect at 3.6MHz.

Above, the calculated Y parameter model including a prototype NT card. This model captures the various loss components of the transformer, mainly magnetising loss, at 3.6MHz. Continue reading Small efficient matching transformer for an EFHW – NEC model at 3.6MHz

## Exploration of feasibility of ground wave comms on 80m VK2OMD-VK1EA

To be useful, S/N needs to exceed 10dB.

Ham mythology has it that ground wave is good for 100km on 80m without much qualification.

## Ambient noise

Ambient noise is very important, it is one factor of the Signal / Noise ratio which sets the limit for communication.

Using ITU-R P.368-9 we find that expected median noise figure in a residential precinct at 3.6MHz is 57dB. Continue reading Exploration of feasibility of ground wave comms on 80m VK2OMD-VK1EA

## nanoVNA – VSWR in terms of 400Ω

I was asked whether the nanoVNA can display VSWR in terms of 400Ω or some other arbitrary impedance.

Some antenna analysers and VNAs support display of results in terms of some specified impedance other than 50Ω, sometimes only a limited fixed set.

The direct answer to the question is “probably no, not directly on the ‘original’ nanoVNA today”, there are many firmware forks and many PC clients for nanovna, and now many significantly hardware versions appearing. Things may change.

However, if you can save a set of impedance measurements, they can be converted to VSWR relative to some other impedance reference.

## An example

At Implementation of G5RV inverted V using high strength aluminium MIG wire – impedance measurements a set of impedance measurements made with a Rigexpert AA600 is presented.

A similar set of measurements could be made with a standalone nanoVNA and saved, but given that it is such a clumsy device and its USB socket has become so unreliable, I will not repeat the measurement but use the data used for the article above.

So we have a table of frequency and (R,X) measured looking into a ‘real’ 450Ω feed line, so we will calculate wrt 450Ω but we could use any reference.

Above is a plot or (R,X) looking into the feed line. Continue reading nanoVNA – VSWR in terms of 400Ω

## Annual inspection of antenna reveals balun box issue

During annual inspection of my HF antenna system, I removed the lid from the common mode choke.

The choke is housed in a (German) Hensel PS electrical terminal box. Hairline cracks were observed radiating out from one of the terminal screws, and adjacent to both screws in the lid. The lid is obviously not strong enough to be secured by two screws in the way done. Continue reading Annual inspection of antenna reveals balun box issue

## Review of W8TEE, AC8GY STL (Radcom Feb 2020)

(Purdum 2020) describes a small transmitting loop (STL) which is a little novel in that it uses an arrangement of four circular conductor loops, two in parallel in series with the other two in parallel.

The article goes on to claim some pretty extraordinary efficiency calculated from radiation resistance for a loop structure that is shown at a height of perhaps 2m above natural ground. Continue reading Review of W8TEE, AC8GY STL (Radcom Feb 2020)

## Stub matching loss can bite you

An article by K2PO in QST Feb 2020 entitled An SWR shifting T illustrates the pitfalls in naive design and implementation of transmission line matching systems. I say naive because the article does not address the matter of loss, yet QST publishes it as an example. Continue reading Stub matching loss can bite you

## A thinking exercise on Jacobi Maximum Power Transfer #4

The article A thinking exercise on Jacobi Maximum Power Transfer #3 discussed Kurokawa’s power reflection coefficient as in indicator of mismatch at a system node.

Above is a demonstration circuit in Simsmith, a linear source with Thevenin equivalent impedance of 50-j5Ω. The equivalent voltage is specified by useZo, which like much of Simsmith is counter intuitive (as you are not actually directly specifying generator impedance):

Vthev and Zthev are chosen so that ‘useZo’ will deliver 1 watt to a circuit impedance that equals the G.Zo. Zthev will be Zo*.


text:
Using Mathematica to see the world in a different light - part I

Posted 5 years ago
11483 Views
|
4 Replies
|
31 Total Likes
|

The international year of light has just drawn to an end. From 4-6 February the closing ceremony took place in Merida, Yukatan.

GeoGraphics[GeoMarker[Entity["City", {"Merida", "Yucatan", "Mexico"}]], GeoRange -> Quantity[3000, "Kilometers"], GeoBackground -> "ReliefMap"]

The international Year of Light was a global initiative of the United Nations to celebrate light and light based technologies. Mathematica's built-in Wikipedia data contains detailed information on the Year of Light; here is the first sentence of the article.

TextSentences[WikipediaData["Year of Light"]][[1]]

The International Year of Light and Light-based Technologies, 2015 (IYL 2015) is a United Nations observance that aims to raise awareness of the achievements of light science and its applications, and its importance to humankind.

I am planning to write three posts to show how the Wolfram Language, its wealth of data, and connected devices can be used to keep the year of light alive at your home. In this first part, I will use a spectrometer, connect it to the Wolfram Language, and try to "see the world in a different light". It turns out that the Wolfram Language will be as important as the hardware for this project.

Light is key to life on earth and sight is a key sense; most of the information our brains process comes from our vision. When the first organism developed primitive vision, that was an enormous evolutionary advantage, allowing them to escape preditors and localise prey. Light plays a crucial role in our modern lifes - most likely the information on this very website was delivered to you using optical fibres and light. Light also allows us to study everything from the smallest particles up to the farthest reaches of the universe.

In the 17th century Sir Isaac Newton introduced the word "spectrum" into optics, referring to the range of colours observed whey light passed through a prism. Today spectrometers are used in many scientific labs to study everything from molecules to the light of stars. In this Community several posts have described the construction of spectrometers and raspberry pi spectrometers. For this blog I will use a commercial, small spectrometer, the C12666MA Micro-Spectrometer, that attaches to an Arduino Uno. Hooking it up to the Arduino is trivial as the pins nicely align.

I will first upload the following code to the Arduino:

// This code is a modified from the original sketch from Peter Jansen
// https://github.com/tricorderproject/arducordermini
// also this code just prints the output to csv output to the terminal.

#define SPEC_GAIN        A0
//#define SPEC_EOS         NA
#define SPEC_ST          A1
#define SPEC_CLK         A2
#define SPEC_VIDEO       A3
#define WHITE_LED        A4
#define LASER_404        A5

#define SPEC_CHANNELS    256
uint16_t data[SPEC_CHANNELS];

void setup() {

//pinMode(SPEC_EOS, INPUT);
pinMode(SPEC_GAIN, OUTPUT);
pinMode(SPEC_ST, OUTPUT);
pinMode(SPEC_CLK, OUTPUT);

pinMode(WHITE_LED, OUTPUT);
pinMode(LASER_404, OUTPUT);
digitalWrite(WHITE_LED, LOW);
digitalWrite(LASER_404, LOW);

//digitalWrite(WHITE_LED, HIGH);
//digitalWrite(LASER_404, HIGH);

digitalWrite(SPEC_GAIN, HIGH);
digitalWrite(SPEC_ST, HIGH);
digitalWrite(SPEC_CLK, HIGH);
digitalWrite(SPEC_GAIN, HIGH); //LOW Gain
//digitalWrite(SPEC_GAIN, LOW); //High Gain

//Serial.begin(9600);
Serial.begin(115200);
}

{
//int delay_time = 35;     // delay per half clock (in microseconds).  This ultimately conrols the integration time.
int delay_time = 1;     // delay per half clock (in microseconds).  This ultimately conrols the integration time.
int idx = 0;
int read_time = 35;      // Amount of time that the analogRead() procedure takes (in microseconds) (different micros will have different times)
int intTime = 5;
int accumulateMode = false;
int i;

// Step 1: start leading clock pulses
for (int i = 0; i < SPEC_CHANNELS; i++) {
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
}

// Step 2: Send start pulse to signal start of integration/light collection
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
digitalWrite(SPEC_ST, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
digitalWrite(SPEC_ST, HIGH);
delayMicroseconds(delay_time);

// Step 3: Integration time -- sample for a period of time determined by the intTime parameter
int blockTime = delay_time * 8;
long int numIntegrationBlocks = ((long)intTime * (long)1000) / (long)blockTime;
for (int i = 0; i < numIntegrationBlocks; i++) {
// Four clocks per pixel
// First block of 2 clocks -- measurement
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);

digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
}

// Step 4: Send start pulse to signal end of integration/light collection
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
digitalWrite(SPEC_ST, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
digitalWrite(SPEC_ST, HIGH);
delayMicroseconds(delay_time);

// Step 5: Read Data 2 (this is the actual read, since the spectrometer has now sampled data)
idx = 0;
for (int i = 0; i < SPEC_CHANNELS; i++) {
// Four clocks per pixel
// First block of 2 clocks -- measurement
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, LOW);

// Analog value is valid on low transition
if (accumulateMode == false) {
} else {
}
idx += 1;

digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);

// Second block of 2 clocks -- idle
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
}

// Step 6: trailing clock pulses
for (int i = 0; i < SPEC_CHANNELS; i++) {
digitalWrite(SPEC_CLK, LOW);
delayMicroseconds(delay_time);
digitalWrite(SPEC_CLK, HIGH);
delayMicroseconds(delay_time);
}
}

void print_data()
{
for (int i = 0; i < SPEC_CHANNELS; i++)
{
Serial.print(data[i]);
Serial.print(',');
}
Serial.print("\n");
}

void loop()
{
//  digitalWrite(LASER_404, HIGH);
//  digitalWrite(LASER_404, LOW);
//  print_data();
//  delay(10);

//  digitalWrite(WHITE_LED, HIGH);
//  digitalWrite(WHITE_LED, LOW);
//  print_data();
//  delay(10);

print_data();
delay(10);
}

The Arduino code is quite long and I first thought to only attach it to the post, but there is a particular line which is important. In the "void readSpectrometer()" part there is the line

int intTime = 5;

which allows you to set the integration time, i.e. the exposure time. That is a very important parameter to get optimal results. With the new features of the Wolfram Language and in particular with the Arduino Device Connection it is possible to adjust the exposure from within the Wolfram Language code. For this post this is not really necessary, so I will save that for another day.

Connecting the Arduino-Spectrometer duo to Mathematica is trivial. Here I connect it to an OSX system.

mySpectrometer = DeviceOpen["Serial", {Quiet[FileNames["tty.usb*", {"/dev"}, Infinity]][[1]], "BaudRate" -> 115200}]

This piece of code

Quiet[FileNames["tty.usb*", {"/dev"}, Infinity]][[1]]

facilitates the connection as it detects the correct device automatically. We can now start collecting data like so:

data = Table[Pause[2]; ToExpression /@ StringSplit[FromCharacterCode[SplitBy[DeviceReadBuffer[mySpectrometer], # == 10 &][[-2]]], ","], {i, 6}];

This data aquisition code actually measures 6 spectra and pauses for 2 seconds between individual measurements. These repeated measurements decrease the noise of the measurements. The measurement procedure is very straight foward. You point the spectrometer at an object and execute the data aquisition. Let's fist look at the spectrum of a fluorescent light bulb.

ListLinePlot[N@Mean[Select[data, Length[#] == 256 &]], PlotRange -> All]

There are several peaks that we will try to understand a little bit later. The x-axis shows 256 bins which represent the different colours, i.e. frequencies. The y-axis shows the count, i.e. the intensity of that frequency band. In oder to be able to interpret the results we first need to calibrate the spectrometer. Here is a link to a table which contains calibration information for various versions of the spectrometer. It turns out that the calibration is performed using a 5th order polynomial; the repective coefficients are given in the calibration table. For my particular spectrometer I obtain:

a0 = 3.170083173*10^2;
b1 = 2.39519817;
b2 = -8.618615345*10^(-4);
b3 = -5.978279712*10^(-6);
b4 = 8.585352787*10^(-9);
b5 = -2.048534811*10^(-12);
wavelength[x_] := a0 + b1  x + b2  x^2 + b3  x^3 + b4*x^4 + b5  x^5

Here is a plot of the calibration curve:

Plot[{wavelength[x]}, {x, 0, 256}]

It transforms the number of the bin to the corresponding wavelength. We can not plot the spectrum with the correct x-axis.

datacalibrated = Transpose@{wavelength /@ Range[256], N@Mean[Select[data, Length[#] == 256 &]]};
ListLinePlot[datacalibrated, PlotRange -> All]

This is much better. This is how a professional spectrum of a fluorescent light bulb looks like:

which is taken from the Wikipedia commons. Here is a list of the peaks and what element they correspond to:

We clearly see the peaks for Mercury, Terbium and Europium. Wolfram|Alpha has a wealth of information about spectral lines and we can use the following line to get it.

WolframAlpha["spectral lines mercury"]

The two dominant lines here are the one at 4046.565 Angstrom and at 4358.335 Angstrom. They correspond to our lines at about 405 nm and 436 nm. Unfortunately, I have failed to extract a list of all relevant spectral lines from Wolfram|Alpha.

Let's try to spice our representation of the spectrum up a bit. It would be nice to have a visual cue as to the colour the different wavelengths correspond to. Mathematica and the Wolfram Language have everything built in to make this task really easy:

ColorData["VisibleSpectrum"][#] & /@ datacalibrated[[All, 1]]

We can now merge this into a band of colours that we can plot with the spectrum.

Graphics[Table[{ColorData["VisibleSpectrum"][datacalibrated[[i, 1]]], Rectangle[{datacalibrated[[i, 1]], 0}, {datacalibrated[[i + 1, 1]],40}]}, {i, 1, Length[datacalibrated] - 1}]]

If we plot this together things start being easier to interpret.

Show[ListLinePlot[datacalibrated, PlotRange -> All], Graphics[Table[{ColorData["VisibleSpectrum"][datacalibrated[[i, 1]]],
Rectangle[{datacalibrated[[i, 1]], 10}, {datacalibrated[[i + 1, 1]], 40}]}, {i, 1, Length[datacalibrated] - 1}]]]

Now there are two more things I want to tweak. First of all there is this off-set of 54; that is a "zero count", i.e. I obtain at least 54 even if there is no signal so I need to subtract that. This number depends on the spectrometer that you have. On a second one I have got that value is different. Also, I would like the lightcurve itself to reflect the colour. The following code achieves that:

Show[ListLinePlot[Evaluate@(Plus[{0., -54.}, #] & /@ datacalibrated),
PlotRange -> {All, {-30, All}}, Joined -> True, Frame -> True,
ColorFunction -> (Blend["VisibleSpectrum", #1* Differences[datacalibrated[[All, 1]][[{1, -1}]]][[1]] +
datacalibrated[[1, 1]]] &), Filling -> Axis, LabelStyle -> Directive[Black, Bold, Medium],
FrameLabel -> {"Wavelength (nm)", "Intensity"}], Graphics[Table[{ColorData["VisibleSpectrum"][datacalibrated[[i, 1]]],
Rectangle[{datacalibrated[[i, 1]], -30}, {datacalibrated[[i + 1, 1]], -10}]}, {i, 1, Length[datacalibrated] - 1}]]]

To identify the spectral lines it is useful to identify maxima of the curve. The following function helps to achieve that:

localMaxPositions =
Compile[{{pts, _Real, 1}},
Module[{result = Table[0, {Length[pts]}], i = 1, ctr = 0}, For[i = 2, i < Length[pts], i++,
If[pts[[i - 1]] < pts[[i]] && pts[[i + 1]] < pts[[i]],result[[++ctr]] = i]];Take[result, ctr]]];

We can now locate the maxima and plot that on the curve:

dplot = ListLinePlot[datacalibrated, PlotRange -> All];
maxs = ListPlot[Select[Nest[#[[localMaxPositions[#[[All, 2]]]]] &, datacalibrated, 1], #[[2]] > 58.33 &], PlotStyle -> Directive[PointSize[0.015], Green]];
Show[{dplot, maxs}]

The following function produces a list of the positons of the maxima:

Select[Nest[#[[localMaxPositions[#[[All, 2]]]]] &, datacalibrated, 1], #[[2]] > 58.33 &][[All, 1]]
(*{366.873694406445,406.4710188519941,436.1860721777116,489.\
5449068238321,544.8439203380104,587.3854766583466,610.\
7844904587837,668.0095922519685,688.8768717778272,707.\
3639428495288}*)

It turns out that there is a packages which contains information about spectral lines, but the elements that we are interested in (Hg, Te, Eu) are not in the database. For elements that are in the database we can compare the measured lines to the ones in the database:

<< ResonanceAbsorptionLines
ElementAbsorptionMap[Na]

What we can do, however, is to generate - or simulate- an approximatio of the spectral lines we measure. There is a discussion on Stackexchange that shows how to plot an emission spectrum. The following lines are taken from that discussion.

spec[wavelength_, width_] := Flatten[Table[{{x, 0, x}, {x, 1, x}}, {x, wavelength - width, wavelength + width, 0.1}], 1];
ListDensityPlot[
spec[#, 1] & /@
Select[Nest[#[[localMaxPositions[#[[All, 2]]]]] &, datacalibrated,
1], #[[2]] > 59.33 &][[All, 1]],
ColorFunction -> ColorData["VisibleSpectrum"],
ColorFunctionScaling -> False, AspectRatio -> .3,
PlotRange -> {{300, 800}},
FrameTicks -> {Automatic, None, None, None},
FrameTicksStyle -> White, Frame -> True, Background -> Black]

Case Studies

Blue sky

I now measure the spectrum of the blue sky. The measurement was taken late in the year and the sky wasn't the "bluest of blues", but we can still make out some interesting features.

mySpectrometer =
DeviceOpen[
"Serial", {Quiet[FileNames["tty.usb*", {"/dev"}, Infinity]][[1]],
"BaudRate" -> 115200}];
Table[Pause[2];
ToExpression /@
StringSplit[
FromCharacterCode[
","], {i, 6}];
a0 = 3.170083173*10^2;
b1 = 2.39519817;
b2 = -8.618615345*10^(-4);
b3 = -5.978279712*10^(-6);
b4 = 8.585352787*10^(-9);
b5 = -2.048534811*10^(-12);
wavelength[x_] := a0 + b1  x + b2  x^2 + b3  x^3 + b4*x^4 + b5  x^5;
datacalibratedsky =
Transpose@{wavelength /@ Range[256], N@Mean[Select[datasky, Length[#] == 256 &]]};
Show[ListLinePlot[Evaluate@(Plus[{0., -54.}, #] & /@ datacalibratedsky),
PlotRange -> {All, {-30, All}}, Joined -> True, Frame -> True, ColorFunction -> (Blend["VisibleSpectrum", #1*Differences[datacalibratedsky[[All, 1]][[{1, -1}]]][[1]] + datacalibratedsky[[1, 1]]] &), Filling -> Axis, LabelStyle -> Directive[Black, Bold, Medium], FrameLabel -> {"Wavelength (nm)", "Intensity"}],
Graphics[Table[{ColorData["VisibleSpectrum"][datacalibratedsky[[i, 1]]], Rectangle[{datacalibratedsky[[i, 1]], -30}, {datacalibratedsky[[i + 1, 1]], -10}]}, {i, 1, Length[datacalibratedsky] - 1}]]]

The spectrum clearly shows the "blue" in the sky. It turns out that the marked double dip corresponds to absorption by water.

Lasers

Using exactly the same code we can now analyse the spectrum of different lasers - red, green and blue. We can clearly see that the red and green lasers have a very narrow peak within the red and green spectral bands, whereas the blue laser has a much broader peak. This is probably related to how blue laser light is generateed in cheap laser pointers.

Here is how the spectrum of a red laser would look like:

ListDensityPlot[
spec[#, 1] & /@
Select[Nest[#[[localMaxPositions[#[[All, 2]]]]] &, datacalibratedlaserr, 1], #[[2]] > 59.33 &][[All, 1]],
ColorFunction -> ColorData["VisibleSpectrum"], ColorFunctionScaling -> False, AspectRatio -> .3,
PlotRange -> {{400, 800}}, FrameTicks -> {Automatic, None, None, None},
FrameTicksStyle -> White, Frame -> True, Background -> Black]

Absorption spectrum - green leaf

Up to now we have maily discussed emission spectra - apart from some features of the blue sky example. Let's now generate an absorption spectrum. In an emission spectrum a material actively emits radiation/light. In an absoprtion spectrum light that passes through the material.

In order to produce a good absorption spectrum we would ideally use a light source that produces a strong continuous spectrum. When I conducted the leaf experiment I only had a very cheap lamp which is usually used by doctors for initial examinations:

I pointed the lamp at a green leaf and measured the light that shone through. The problem was that the lamp did produce a very poor spectrum heavily biased towards the red frequency range. So in this case study I will first measure the emission spectrum of the lamp and then normalise the absorption spectrum of the leaf. Let's start with the lamp.

mySpectrometer =
DeviceOpen[
"Serial", {Quiet[FileNames["tty.usb*", {"/dev"}, Infinity]][[1]],
"BaudRate" -> 115200}];
datalamp =
Table[Pause[2];
ToExpression /@
StringSplit[
FromCharacterCode[
","], {i, 6}];
a0 = 3.170083173*10^2;
b1 = 2.39519817;
b2 = -8.618615345*10^(-4);
b3 = -5.978279712*10^(-6);
b4 = 8.585352787*10^(-9);
b5 = -2.048534811*10^(-12);
wavelength[x_] := a0 + b1  x + b2  x^2 + b3  x^3 + b4*x^4 + b5  x^5;
datacalibratedlamp =
Transpose@{wavelength /@ Range[256], N@Mean[Select[datalamp, Length[#] == 256 &]]};
Show[ListLinePlot[Evaluate@(Plus[{0., -54.}, #] & /@ datacalibratedlamp),
PlotRange -> {All, {-30, All}}, Joined -> True, Frame -> True, ColorFunction -> (Blend["VisibleSpectrum", #1*Differences[datacalibratedlamp[[All, 1]][[{1, -1}]]][[1]] + datacalibratedlamp[[1, 1]]] &), Filling -> Axis,
LabelStyle -> Directive[Black, Bold, Medium],
FrameLabel -> {"Wavelength (nm)", "Intensity"}],
Graphics[Table[{ColorData["VisibleSpectrum"][datacalibratedlamp[[i, 1]]],
Rectangle[{datacalibratedlamp[[i,1]], -30}, {datacalibratedlamp[[i + 1, 1]], -10}]}, {i, 1,Length[datacalibratedlamp] - 1}]]]

The spectrum is clearly biased to the red freqencies. In the blue/ultraviolet range the light source produces to little output that it will be impossible to determine a reasonable absorption spectrum there. Let's go on to measure the absorption of the leaf.

absorptionleaf =
Transpose[{datacalibratedleaf[[All, 1]], (datacalibratedleaf[[All, 2]]/datacalibratedlamp[[All, 2]])}];
Show[ListLinePlot[Evaluate@(Plus[{0., 0.}, #] & /@ absorptionleaf), PlotRange -> {All, {-0.3, 1}}, Joined -> True, Frame -> True,
ColorFunction -> (Blend["VisibleSpectrum", #1*Differences[absorptionleaf[[All, 1]][[{1, -1}]]][[1]] + absorptionleaf[[1, 1]]] &), Filling -> Axis,
LabelStyle -> Directive[Black, Bold, Medium], FrameLabel -> {"Wavelength (nm)", "Intensity"}], Graphics[Table[{ColorData["VisibleSpectrum"][absorptionleaf[[i, 1]]],Rectangle[{absorptionleaf[[i, 1]], -0.3}, {absorptionleaf[[i + 1,1]], -.1}]}, {i, 1, Length[absorptionleaf] - 1}]]]

The absorption in the red frequency range stems from chlorophyll; chlorophyll also absorbs in the 400-450nm range which is hard to se here, because of our poor light source. The dip at around 500nm is carotenoids. You can compare that to the spectrum of leaves on this website.

"Black body radiation" - an old fashioned light bulb

The final example will be of an old fashioned, small light bulb.

I will use a bench supply to slowly increase the voltage and the current. The filament will go from a red glowing colour to brighter "whiter" colour, but we will see that even for the highest voltage of 12V the spectrum will still be quite different from "white", i.e. uniform. The code is just the same as above. I save all plots for voltages 1V to 12V with increments of 1V in one variable:

specall = {spec1V, spec2V, spec3V, spec4V, spec5V, spec6V, spec7V, spec8V, spec9V, spec10V, spec11V, spec12V}

This can be easily plotted like so:

Grid[Partition[specall, 4]]

We can also animate this:

ListAnimate[specall, AnimationRunTime -> 10]

Note that from 9V onwards the spectrum has max-ed out. When I measured the spectra I took note of the voltage and the corresponding current:

voltageamp = {{1, 0.04}, {2, 0.06}, {3, 0.07}, {4, 0.08}, {5, 0.1}, {6, 0.11}, {7, 0.12}, {8, 0.13}, {9, 0.14}, {10, 0.14}, {11, 0.15}, {12, 0.16}}

I can now represent that with the regression line.

Show[ListPlot[voltageamp, AxesLabel -> {"Volts", "Amps"}, LabelStyle -> Directive[Bold, Medium]],
Plot[Evaluate@Fit[voltageamp, {1, x}, x], {x, 0, 12}, PlotStyle -> Red]]

Note the nearly linear increase in current as the voltage increases. Finally, I took photos as the voltage increased.

We now add a little bit of motion to the voltage amp graph:

figvoltsamp =
Evaluate@Table[Show[ListPlot[voltageamp, AxesLabel -> {"Volts", "Amps"}, LabelStyle -> Directive[Bold, Medium],
Epilog -> {PointSize[Large], Green, Point[voltageamp[[i]]]}], Plot[Evaluate@Fit[voltageamp, {1, x}, x], {x, 0, 12}, PlotStyle -> Red]], {i, 1, 12}];

To finish everything off, we can now animate this:

ListAnimate[GraphicsRow[#, ImageSize -> Full] & /@ Transpose[{specall, bulb, figvoltsamp}]]

If we were to assume that this is the radiation of a black body, we could use Planck's law:

FormulaData[{"PlanckRadiationLaw", "Wavelength"}]

We could plot this as a function of temperature and wavelength

equation =
"\[Lambda]" -> Quantity[l, "Nanometers"]}];
Plot3D[Quantity[1.191042*^29, ("Pascals")/("Seconds")]/(-1.16.255 l^5 + 2.71828^(1.438*^7/(l t)) l^5), {l, 300, 808}, {t, 1000, 5000}, PlotRange -> All,
AxesLabel -> {"wavelength", "temperature", "luminosity"}, LabelStyle -> Directive[Bold, Medium], ImageSize -> Large]

We can also write a little loop to calculate for which wavelength the maximum for different temperatures is reached.

results = {}; Monitor[
"Wavelength"}, {"T" -> Quantity[t, "Kelvins"], "\[Lambda]" -> Quantity[l, "Nanometers"]}];
root = FindRoot[D[equation2[[2, 2]], l] == 0, {l, 700}, MaxIterations -> Infinity];
AppendTo[results, {t, root}], {t, 2800, 3400, 100}], t]

This gives the following table:

Grid[Join[{{"Temperature K", "\[Lambda] at max"}}, Transpose[{results[[All, 1]], results[[All, 2, -1, 2]]}]], Frame -> All]

We can plot the relationship between maxium and temperature.

ListLinePlot[Transpose[{results[[All, 1]], results[[All, 2, -1, 2]]}],Mesh -> Full, MeshStyle -> Red,
AxesLabel -> {"Temperature", "wavelength at max"}, LabelStyle -> Directive[Bold, Medium], ImageSize -> Large]

This allows us in principle to estimate the temperature of the filament. Note, that for 8V the maximum is at about 750nm. The graph shows that this corresponds to a temperature higher than 3400K, which is too high for such a small light bulb. So there is still quite some room for improvement... The principle of temperature measurement, by the colour of the sample is realistic.

There are many further projects one could think of. I suppose that with a decent telescope it should be possible to analyse the light of stars for example. Also, the spectrometer works nicely on a Raspberry Pi. It should be quite straight forward to take measurements of say the sky over the day and see how the dominant colours change.

Cheers,

Marco

4 Replies
Sort By:
Posted 5 years ago
 Awesome post, Marco! I'm looking forward to reading your other two posts to "see the world in a different light".
Posted 5 years ago
 Dear @Vitaliy Kaurov ,thank you very much for your kind and encouraging words. I am glad that you liked the post. I was, in fact, aware that one can use Wolfram|Alpha and then click on the more button in the subpod. My problem is that the data I managed to download appears to be incomplete. If we use Wolfram|Alpha to display the spectrum of say mercury it looks like this:There are spectral lines all over the visual range from 4000-7500 Angstroms. If I go to the pod content I only get the first lines. I can click on "More" and I get some more lines, but no matter how often I click on "More" I only get 21 lines: lines=WolframAlpha["Spectral lines mercury", {{"Lines:AtomicSpectrumData", 1}, "ComputableData"}, PodStates -> {"Lines:AtomicSpectrumData__More"}] So in a table this looks like this: lines[[1, 1]] // TableForm The last lines I download are at about 4700 Angstroms. That means that a large part of the visual spectrum is not covered in the list I download. I have tried various modifications of the Wolfram|Alpha request, but nothing downloaded all the visible lines. This is not a problem for elements that only have "few" spectral lines, but is more of a problem for those, like mercury, terbium etc that have lots of them. I am not terribly familiar with Wolfram|Alpha syntax and do not know whether there is a way to extract a list of all lines. I am aware that many of the lines are quite weak. Wolfram Alpha also has information on the intensity of the lines: linesplus = WolframAlpha["spectrum mercury", {{"Lines:AtomicSpectrumData", 1}, "ComputableData"}, PodStates -> {"Lines:AtomicSpectrumData__Show line \characteristics"}]; linesplus[[1, 1]] // TableForm The transition probability tells us which lines would be strong, i.e. strong enough to be seen with a cheap spectrometer; it appears that only the ones with transition probabilities of about $10^7$ are really easy to detect. I have also noted that the same requests I used above for mercury do not work equally well for all elements, say terbium. The objective would be to input a spectrum, like the ones I have measured in the post. They I would like to compare that to a database of known spectra and program something that attempts to guess at what is contained in the sample. I have tried importing data from some spectral databases, like this one or this one or this one , but nothing really worked well enough. It would be nice to be able to detect these components from the spectrum. In a project we are doing here, we try to send the spectral data to mobile phones. If we could use CloudDeploy on these platforms, perhaps even on the Apple Watch, that would be nice.We have seen many fantastic posts here where devices were hooked up to Mathematica and where data was read. Most of the time Mathematica is "only" used to represent data, stream the data online, read it from other deices etc. I wonder whether: (i) the Wolfram Language's abilities could be used even more to improve the data from poor devices, e.g. perform a sort of deconvolution on the spectra I got or do some tricks to mathematically improve the measurements, and (ii) I am looking for examples of how the Wolfram Language adds substantially to the analysis of the data. Being able to decompose the spectrum into the components would be a nice case study. As I said, I will try to post at least two further examples and I hope that I'll be able to show how much the Wolfram Language can improve the data from poor sensors and perhaps more...I'd love to hear your thoughts on this.Best wishes,Marco
 @Marco Thiel, amazing post as usual, and we excitingly adding it to Staff Picks! At some point you said: I have failed to extract a list of all relevant spectral lines from Wolfram|Alpha Did you mean the following code? mercuryLINES = WolframAlpha["spectral lines mercury", {{"Lines:AtomicSpectrumData", 1}, "ComputableData"}]; TableForm[Rest[#[[1, 1]]], TableHeadings -> {None, First[#[[1, 1]]]}] &@mercuryLINES You do not need to know the syntax of WolframAlpha function, just access it interactively (as explained here: Data Formats in Wolfram|Alpha):


text:
# Multiple Equilibria in a Land–Atmosphere Coupled System

• Corresponding author: Jianping HUANG, hjp@lzu.edu.cn
• Funds:

Supported by the National Science Foundation of China (41521004 and 41705047), Strategic Priority Research Program of Chinese Academy of Sciences (XDA2006010301), Foundation of Key Laboratory for Semi-Arid Climate Change of the Ministry of Education in Lanzhou University from the Fundamental Research Funds for the Central Universities (lzujbky-2017-bt04), and China 111 Project (B 13045)

• doi: 10.1007/s13351-018-8012-y
• Many low-order modeling studies indicate that there may be multiple equilibria in the atmosphere induced by thermal and topographic forcings. However, most work uses uncoupled atmospheric model and just focuses on the multiple equilibria with distinct wave amplitude, i.e., the high- and low-index equilibria. Here, a low-order coupled land–atmosphere model is used to study the multiple equilibria with both distinct wave phase and wave amplitude. The model combines a two-layer quasi-geostrophic channel model and an energy balance model. Highly truncated spectral expansions are used and the results show that there may be two stable equilibria with distinct wave phase relative to the topography: one (the other) has a lower layer streamfunction that is nearly in (out of) phase with the topography, i.e., the lower layer ridges (troughs) are over the mountains, called ridge-type (trough-type) equilibria. The wave phase of equilibrium state depends on the direction of lower layer zonal wind and horizontal scale of the topography. The multiple wave phase equilibria associated with ridge- and trough-types originate from the orographic instability of the Hadley circulation, which is a pitch-fork bifurcation. Compared with the uncoupled model, the land–atmosphere coupled system produces more stable atmospheric flow and more ridge-type equilibrium states, particularly, these effects are primarily attributed to the longwave radiation fluxes. The upper layer streamfunctions of both ridge- and trough-type equilibria are also characterized by either a high- or low-index flow pattern. However, the multiple wave phase equilibria associated with ridge- and trough-types are more prominent than multiple wave amplitude equilibria associated with high- and low-index types in this study.
• Fig. 1.  The second one (left panels) and third one (right panels) of the three equilibrium states for m = 3.7 at Cg = 50 W m–2. They belong to “High 2” and “High 1” equilibria, respectively. The streamfunction fields of the (a, e) upper and (b, f) lower layers, respectively. The temperature fields of (c, g) the atmosphere and (d, h) the land, respectively. The contour intervals are (a, e) 2.0 × 107 m2 s–1, (b) 2.0 × 105 m2 s–1, (f) 3.0 × 105 m2 s–1, and (c, d, g, h) 10 K. The background dotted lines show the topographic heights in the model, with negative regions shaded.

Fig. 2.  As in Fig. 1, but for the third one of the three equilibrium states for m = 3.7 at Cg = 55 W m–2 (left panels) and for the third one of the five equilibrium states for m = 6 at Cg = 30 W m–2 (right panels). They belong to “Low 1” and “Low 2” equilibria, respectively. The “Low 2” equilibria has nondimensional solutions with $({\psi _1}, \;{\psi _2}, \;{\psi _3}, \;{\theta _1}, \;{\theta _2}, \;{\theta _3}, \;{T_{{\rm g}, 1}}, \;{T_{{\rm g}, 2}}, \;{T_{{\rm g}, 3}})$ = (0.0216, –0.0079, –0.0051, 0.0311, –0.0098, –0.0043, 0.0771, –0.0175, –0.0077). The contour intervals are (a) 2.0 × 107 m2 s–1, (e) 1.0 × 107 m2 s–1, (b, f) 1.0 × 106 m2 s–1, (c, d) 10 K, and (g, h) 4 K.

Fig. 3.  The equilibrium bifurcation associated with the change in meridional differential solar heating parameter Cg for m = 3.7 (left panels) and m = 6 (right panels), respectively. The ordinate shows the nondimensional equilibrium values of (a, d) the zonal component $\psi _1^1$ and the wave components (b, e) $\psi _2^1$ and (c, f) $\psi _3^1$, respectively. Different branches of equilibrium solutions have different colors. The crosses denote unstable equilibria, the circles denote stable high-index equilibria, and the asterisks denote stable low-index equilibria.

Fig. 5.  (a) Stability curves of the Hadley circulation in the coupled land–atmosphere model (Case 1). The blue solid lines enclose the region of orographic instability. The red solid (black dashed) lines and the top x-axis and the right y-axis enclose the region of baroclinic instability in the presence (absence) of topography. Comparison of the regions of (b) orographic instability, (c) baroclinic instability in the presence of topography, and (d) baroclinic instability in the absence of topography for the four experiments (Cases 1–4).

Fig. 4.  A traveling wave solution in the absence of topography for m = 3.7 at Cg = 50 W m–2, with $({\psi _1}, \;{\psi _2}, \;{\psi _3}, \;{\theta _1}, \;{\theta _2}, \;{\theta _3}, \;{T_{{\rm g}, 1}}, \;{T_{{\rm g}, 2}}, \;{T_{{\rm g}, 3}})$= (0.0534, 0.0214, 0.0016, 0.0534, 0.0125, –0.0008, 0.1311, 0.0014, 0.0048) at this moment. (a, b) The streamfunction fields of the upper and lower layers, respectively. The contour intervals are (a) 2.0 × 107 m2 s–1 and (b) 2.0 × 106 m2 s–1. (c) Temporal evolution for the nondimensional equilibrium values of the zonal component $\psi _1^1$ (solid line) and $\psi _1^3$ (dashed line). (d) Temporal evolution for the nondimensional equilibrium values of the wave component $\psi _3^1$ (solid line) and $\psi _3^3$ (dashed line).

Fig. 6.  As in Fig. 3, but for Case 3 (without heat flux) with (a) m = 3.7 and (b) m = 6, (c) for Case 2 (without coupling) with m = 3.7, and (d) for Case 4 (without longwave radiation) with m = 3.7. Each ordinate shows the nondimensional equilibrium solution of the wave component $\psi _3^1$.

Fig. 7.  As in Fig. 1, but for the only stable equilibrium state for m = 3.7 at Q = 50 W m–2 in Case 2 (left panels) and for the third one of the three equilibrium states for m = 3.7 at Cg = 50 W m–2 in Case 4 (right panels). The former has nondimensional solutions with $({\psi _1}, \;{\psi _2}, \;{\psi _3}, \;{\theta _1}, \;{\theta _2}, \;{\theta _3})$ = (0.0656, –0.0679, 0.0259, 0.0613, –0.0616, 0.0253). The latter has nondimensional solutions with $({\psi _1}, \;{\psi _2}, \;{\psi _3}, \;{\theta _1}, \;{\theta _2}, \;{\theta _3}, \;{T_{{\rm g}, 1}}, \;{T_{{\rm g}, 2}}, \;{T_{{\rm g}, 3}})$ = (0.0662, –0.0863, 0.0278, 0.0615, –0.0776, 0.0272, 0.1762, –0.1551, 0.0544). They both belong to “Low 1” equilibria. The contour intervals are (a, e) 4.0 × 107 m2 s–1, (b, f) 2.0 × 106 m2 s–1, and (c, g, h) 10 K. Note that there is no land temperature field in Case 2 (left panel).

Fig. 8.  As in Fig. 1, but for the second one (left panels) and third one (right panels) of the three equilibrium states for m = 3.7 at Cg = 50 W m–2 in Case 3. They have nondimensional solutions with $({\psi _1}, \;{\psi _2}, \;{\psi _3}, \;{\theta _1}, \;{\theta _2}, \;{\theta _3}, \;{T_{{\rm g}, 1}}, \;{T_{{\rm g}, 2}}, \;{T_{{\rm g}, 3}})$ = (0.0623, –0.0002, –0.0020, 0.0626, –0.0004, –0.0019, 0.1596, –0.0004, –0.0021) and (0.0628, –0.0006, 0.0052, 0.0621, –0.0003, 0.0052, 0.1591, –0.0003, 0.0056), respectively. They belong to “High 2” and “High 1” equilibria, respectively. The contour intervals are (a, e) 2.0 × 107 m2 s–1, (b, f) 1.0 × 105 m2 s–1, and (c, d, g, h) 10 K.

Fig. 9.  The heating fields of the “High 2” (left panels) and “High 1” (right panels) equilibrium states shown in Fig. 1, respectively. (a, e) The shortwave radiation, (b, f) the longwave radiation, (c, g) the heat flux, and (d, h) the net diabatic heating absorbed by the atmosphere. All of the contour intervals are 10 W m–2. The background dotted lines show the topographic heights in the model, with negative regions shaded.

Fig. 10.  (a)–(d) The net diabatic heating absorbed by the atmosphere for the equilibrium states shown in Figs. 7 and 8, respectively. The contour intervals are (a, b) 30 W m–2 and (c, d) 5 W m–2. The background dotted lines show the topographic heights in the model, with negative regions shaded.

Fig. 11.  As in Fig. 3, but for dimensional variables for m = 3.7. Only stable equilibrium states are shown here. The blue (red) branch represents trough-type (ridge-type) equilibria, and the black branch represents the Hadley equilibria.

Fig. 12.  As in Fig. 11, but for m = 6.

Fig. 13.  Phase diagrams of dimensional variables for m = 3.7 (left panels) and m = 6 (right panels), respectively. Each ordinate shows the variable AH, and the abscissa gives (a, e) ΔTa, (b, f) ΔTg, (c, g) ATa, and (d, h) ATg. The meaning of colors and symbols are same as that in Fig. 3.

###### 通讯作者: 陈斌, bchen63@163.com
• 1.

沈阳化工大学材料科学与工程学院 沈阳 110142

## Multiple Equilibria in a Land–Atmosphere Coupled System

###### Corresponding author: Jianping HUANG, hjp@lzu.edu.cn;
• Key Laboratory for Semi-Arid Climate Change of the Ministry of Education, College of Atmospheric Sciences, Lanzhou University, Lanzhou 730000
Funds: Supported by the National Science Foundation of China (41521004 and 41705047), Strategic Priority Research Program of Chinese Academy of Sciences (XDA2006010301), Foundation of Key Laboratory for Semi-Arid Climate Change of the Ministry of Education in Lanzhou University from the Fundamental Research Funds for the Central Universities (lzujbky-2017-bt04), and China 111 Project (B 13045)

Abstract: Many low-order modeling studies indicate that there may be multiple equilibria in the atmosphere induced by thermal and topographic forcings. However, most work uses uncoupled atmospheric model and just focuses on the multiple equilibria with distinct wave amplitude, i.e., the high- and low-index equilibria. Here, a low-order coupled land–atmosphere model is used to study the multiple equilibria with both distinct wave phase and wave amplitude. The model combines a two-layer quasi-geostrophic channel model and an energy balance model. Highly truncated spectral expansions are used and the results show that there may be two stable equilibria with distinct wave phase relative to the topography: one (the other) has a lower layer streamfunction that is nearly in (out of) phase with the topography, i.e., the lower layer ridges (troughs) are over the mountains, called ridge-type (trough-type) equilibria. The wave phase of equilibrium state depends on the direction of lower layer zonal wind and horizontal scale of the topography. The multiple wave phase equilibria associated with ridge- and trough-types originate from the orographic instability of the Hadley circulation, which is a pitch-fork bifurcation. Compared with the uncoupled model, the land–atmosphere coupled system produces more stable atmospheric flow and more ridge-type equilibrium states, particularly, these effects are primarily attributed to the longwave radiation fluxes. The upper layer streamfunctions of both ridge- and trough-type equilibria are also characterized by either a high- or low-index flow pattern. However, the multiple wave phase equilibria associated with ridge- and trough-types are more prominent than multiple wave amplitude equilibria associated with high- and low-index types in this study.

Reference (51)

/


text:
Firefox failing to find plugin when path too long >= 130 chars

RESOLVED INCOMPLETE

Status

()

Core
Plug-ins
--
major
RESOLVED INCOMPLETE
7 years ago
23 days ago

1.9.2 Branch
---
x86
Windows XP
regression
Points:
---

Details

(Reporter)

Description

7 years ago
User-Agent:       Mozilla/5.0 (Windows; U; Windows NT 5.1; en-US) AppleWebKit/533.4 (KHTML, like Gecko) Chrome/5.0.375.99 Safari/533.4
Build Identifier: Mozilla/5.0 (Windows; U; Windows NT 5.1; en-US; rv:1.9.2.6) Gecko/20100625 Firefox/3.6.6

Firefox is not finding all of our plugins.
The same plugins are found in recent, older versions of Firefox.

Other versions tested:
3.5.5       - works
3.5.10      - works
3.6.4       - fails
4.0 beta 1  - fails

Symptoms:
- Some plugins are missing from about:plugins.
- When adding an embed tag to a page using JavaScript the browser displays the
missing plugin UI

Hypothesis:
After comparing the windows registry settings of a plugin that was found vs
one that wasn't it appears the browser may be failing to find the plugin
becuase of a path length limitation.
The steps to reproduce in this report seem to validate this hypothesis.

Reproducible: Always

Steps to Reproduce:
1.Edit windows registry for a given plugin to use a path >= 130 characters (us ascii) in length.
Registry key\value used in these steps is:
[HKEY_CURRENT_USER\Software\MozillaPlugins\@your_org/your_plugin_desc;version=a.b.c.d]\Path="absolute_path\\to\\your\\npYourPlugin.dll"

2. Place the plugin at the path you specified in the registry in #1 above (>= 130 chars).

3. Start firefox and observe that the plugin does NOT appear in about:plugins

4. Close firefox, re-edit the windows registry using a path with < 130 characters.

5. Restart firefox and observe that plugin will now appear in about:plugins.

6. Yay - you reproduced the problem - get coffee now :)
Actual Results:
The browser will find plugins, as observed in about:plugins when the path to the plugin DLL is < 130 characters and will not find them when it is >= 130.

Expected Results:
Firefox should be able to support the max path length for the filesystem from which it searches for and loads plugins.
In past versions (see summary in this report) it supported longer paths.
(Reporter)

Updated

7 years ago
Version: unspecified → 3.6 Branch

Updated

7 years ago
Component: Extension Compatibility → Plug-ins
Product: Firefox → Core
QA Contact: extension.compatibility → plugins
Version: 3.6 Branch → 1.9.2 Branch

Updated

7 years ago
Keywords: regression, regressionwindow-wanted
(Reporter)

Comment 1

7 years ago
re - regression window:
Problem started sometime between 3.5.10 and 3.6.4.

Comment 2

7 years ago
Similar, but unrelated to this is a limit on the maximum number of characters in the registry key to less than the 'real' maximum for a single API call is 32 * (255 + 1) characters [32 levels deep w/ 255 key length limit + \ separator].

Dynamic length calculation should be used alongside malloc, not _MAX_PATH which is not correct in Windows, especially when Unicode is involved.

Problem #1:
_MAX_PATH is 260 in Windows, whereas we are hitting a 130 character limit...

This is due to the fact that NS_ARRAY_LENGTH is returning the number of elements (characters) in the array.  RegQueryValueExW expects that the length passed into it be the number of bytes in the buffer... even though the output is in wide characters.

Changing:

DWORD pathlen = NS_ARRAY_LENGTH(path);
to
DWORD pathlen = NS_ARRAY_LENGTH(path) * sizeof(path[0]);

Should at least bring the length to 260.

Problem #2: Assuming _MAX_PATH in Windows is "really" the maximum path length.  This may (hopefully) hold for many function in the C API that Windows exports, especially ones that don't accept an input length on the buffer to be written to (ex: tmpnam).
The retrieval of the path from the registry should query the length of the value, then use that to get a properly sized buffer to store the full path in.
Querying the length of the path would be done through passing a NULL data buffer in, but passing the address of the data buffer length variable.

Comment 3

7 years ago
This regressed on MC within:
http://hg.mozilla.org/mozilla-central/pushloghtml?fromchange=44c392db6672&tochange=96e8d529b2d3
This includes an Electrolysis Branch Merge to MC.
Hence i tried to find a Range on that Branch, but even going back to the first working Windows Build (Built from http://hg.mozilla.org/projects/electrolysis/rev/ef5eac530f66) the Issue is given.
Status: UNCONFIRMED → NEW
Ever confirmed: true
Keywords: regressionwindow-wanted

Comment 4

7 years ago
(In reply to comment #2)
> DWORD pathlen = NS_ARRAY_LENGTH(path) * sizeof(path[0]);
It's redundant because the preprocessor will expand this to
"sizeof(path) / sizeof(path[0]) * sizeof(path[0])".

We don't have to consider about "real" maximum path length here because we already depends on MAX_PATH in many places.
http://mxr.mozilla.org/mozilla-central/search?string=MAX_PATH
Fixing this will be too big to land on the branch for the regression.

Comment 5

6 years ago
I understand about not considering about "real" maximum path length. Also, in depoyment to end user machine, it should usually not be much of a problem. Except when the plugin is deployed for an individual non-admin user instead of the whole machine. In this case, the path can get closer to the limit. For example Google update can have a path like "C:\Documents and Settings\gmzwinge\Local Settings\Application Data\Google\Update\1.2.183.39\npGoogleOneClick8.dll" which is 114 characters. Also, while developing plugin, that path length can be more easily reached due to the multiple hierarchy of folders usually involved in development. Would it be possible to make a change similar to the one proposed by Thomas above? The change goes in mozilla-1.9.\modules\plugin\base\src\nsPluginDirServiceProvider.cpp, in the function nsPluginDirServiceProvider::GetPLIDDirectoriesWithHKEY(). The change could be from:
WCHAR path[_MAX_PATH];
DWORD pathlen = NS_ARRAY_LENGTH(path);
to:
WCHAR path[_MAX_PATH];
DWORD pathlen = sizeof(path);
This would double the length of the path. And hopefully should not be too big of a change. Thanks.

Comment 6

5 years ago
Hi Guys, just to update you, finally I found a way to solve it
Do you have error messages?, like :

* Path too long
* Error cannot delete file: cannot read from source file or disk
* Cannot delete file: Access is denied
* There has been a sharing violation.
* Cannot delete file or folder The file name you specified is not valid or too long. Specify a different file name.

For that I tried it with: http://longpathtool.com/

Comment 7

23 days ago
Resolving old bugs which are likely not relevant any more, since NPAPI plugins are deprecated.
Status: NEW → RESOLVED
Last Resolved: 23 days ago
Resolution: --- → INCOMPLETE
You need to log in before you can comment on or make changes to this bug.