Coding

## Julia 101 – 0.1. Least-Square fitting

``````year = Float64[0]

for i in 1:6
push!(year, i)
end

year = year .+ 2013``````
```7-element Array{Float64,1}:
2013.0
2014.0
2015.0
2016.0
2017.0
2018.0
2019.0```
``````year_x = collect(Float64, 2013:1:2020);
h = [2, 2, 2, 3, 4, 6, 8, 12];``````
``````using Plots
scatter(year_x, h,
xlab = "year", ylab = "h-index",
xlim = [2012, 2023], ylim = [0,15], framestyle = :box)

scatter!(year_x[3:end], h[3:end])``````
``````using LsqFit

@. model(x, p) = p[1] + x*p[2] + x^2.0*p[3]
fit = curve_fit(model, year_x, h, [2.0, 0.5, -3.0])``````
`LsqFit.LsqFitResult{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Int64,1}}([1.1832864618761656e6, -1174.946148708822, 0.2916665972735931], [0.20833284547552466, -0.19642864260822535, -0.017856936203315854, -0.2559520348440856, 0.08928606053814292, 0.017857350641861558, 0.5297618352342397, -0.3750004852190614], [1.0000000000082452 2012.9999999932165 4.05216899997484e6; 1.0000000000082452 2014.0000000073403 4.05619599998367e6; … ; 1.0000000000082452 2018.9999999797847 4.076360999975335e6; 1.0000000000082452 2019.9999999939084 4.08039999999393e6], true, Int64[])`
``````T = typeof(fit)
for (name, typ) in zip(fieldnames(T), T.types)
println("type of the fieldname typ")
end``````
```type of the fieldname param is Array{Float64,1}
type of the fieldname resid is Array{Float64,1}
type of the fieldname jacobian is Array{Float64,2}
type of the fieldname converged is Bool
type of the fieldname wt is Array{Int64,1}
```
``?fit``
```search: fit LsqFit filter filter! curve_fit first firstindex isfinite popfirst!

```

No documentation found.

`fit` is of type `LsqFit.LsqFitResult{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Int64,1}}`.

# Summary

``struct LsqFit.LsqFitResult{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Int64,1}} &lt;: Any``

# Fields

``````param     :: Array{Float64,1}
resid     :: Array{Float64,1}
jacobian  :: Array{Float64,2}
converged :: Bool
wt        :: Array{Int64,1}``````
``fitted_p = fit.param``
```3-element Array{Float64,1}:
1.1832864618761656e6
-1174.946148708822
0.2916665972735931```
``````new_x = collect(Float64, 2013:1:2025)
new_y = fitted_p[1] .+ new_x.*fitted_p[2] .+ (new_x.^2).*fitted_p[3]
plot!(new_x, new_y, ylim = [0,30], seriestype =:line)``````
``fit.resid``
```8-element Array{Float64,1}:
0.20833284547552466
-0.19642864260822535
-0.017856936203315854
-0.2559520348440856
0.08928606053814292
0.017857350641861558
0.5297618352342397
-0.3750004852190614```
``fit.jacobian``
```8×3 Array{Float64,2}:
1.0  2013.0  4.05217e6
1.0  2014.0  4.0562e6
1.0  2015.0  4.06022e6
1.0  2016.0  4.06426e6
1.0  2017.0  4.06829e6
1.0  2018.0  4.07232e6
1.0  2019.0  4.07636e6
1.0  2020.0  4.0804e6```
``fit.wt``
`0-element Array{Int64,1}`
``fit.converged``
`true`

## 利用BASH脚本运行IDL程序、创建子目录或删除文件

Here are the things I want to do:

1.  I have a data folder containing many sources, like Arp220, NGC981, etc.. Each of them possesses a folder under the whole data folder (./source/Arp220, ./source/NGC891, …). I want to create a subfolder named “plot” under each sources folder, e.g., ./source/Arp220/plot.
2. And then I want to copy the IDL routine to each “plot” folder I just created and run all the PRO routines at one time.

Here is a bash script being able to do such things:

```#!/bin/bash
for i in \$(ls -d */); do mkdir -p \$i/'plot'; done
for i in \$(ls -d */); do cp ./r_name.pro \$i/'plot'/; done
for i in \$(ls -d */); do cd \$PWD/\$i/'plot';idl -e ".r ./r_name.pro"; cd ../.. ; done```

The first line means creating a “plot” folder under all the subfolders under you present folder. And the next line is to copy the IDL routine (r_name.pro) to each “plot” folder you created. The last line runs all the IDL routines.

`find ./* -type f|grep X|while read file; do rm "\$file"; done`

`find ./* -type f -name "*X*" -exec rm {} \;`

## Obfuscated ssh Client under Linux

SSH is a powerful tool not only for astronomers. We often use it to download observation data from the observatory’s server or run some complicated programs on remote (overseas) servers. But for the reason as we all know, the SSH tunnels are always being obstructed. So we need some technique to get over such obstruction. Therefore, the obfuscated SSH turns out to be a sufficient method.

I will describe how to use obfuscated SSH_TUNNEL under Linux as below:
(Original Code: Here)

INSTALLATION:

```wget http://aenes.googlecode.com/files/brl-obfuscated-openssh.zip
unzip brl-obfuscated-openssh.zip
cd ./brl-obfuscated*
./configure  --prefix=/usr/local/newssh --sysconfdir=/etc/newssh
make
make install```

Then you need to configure the ssh server following the guide here.

USAGE:
Then use this to run the ssh client:

`/usr/local/newssh/bin/ssh -N -v -Z ObfuscateKeyWord -p ***** username@hostname -D 127.0.0.1:7070`

There are 4 parameters you need to specify: “ObfuscateKeyWord” is the obfuscated keyword you set in the configuration file (sshd_config) of the server, “*****” is the ssh server port number, “username” and “hostname” are as the name described.

## Least-Squares Fitting of data containing errors

(This is a note personal only! )

There is a powerful tool dealing with data containing both “X” and “Y” errors to make a lease-square fitting in IDL, that is mpfit. This routine allows you to use a self-defined function to do the fitting. Also, you can set up some conditions such as the parameters’ limit or fixed parameters etc. This is fabulous! So beginning with the actual astronomical data have errors both in “X” and “Y”, you need to use a function called “linfitex” (you can see more info from the mpfit website) to both define the function you want to fit and calculate the least-square needed in the fitting. Here is an example of the function:

```FUCNTION myfunct, p, x=x, y=y \\$
, sigma_x=sigma_x, sigma_y=sigma_y, _EXTRA=extra
a = p[0]
b = p[1]
f = a*alog10(X)+b
;#########calculate the derivative of the function#########
;##put into resid as (y-f)/sqrt(sigma_y^2+(y'*sigma_x)^2)##
resid = (y - f)/sqrt(sigma_y^2 + (a/(x*alog(10))*sigma_x)^2)
RETURN, resid```

In the code above I just defined a function as

Then you can read your data. and when you want to fit the data, just using mpfit as:

```p      = [1,1,1]
p[1]   = 1
p[2]   = 2
p[3]   = 3 ;initial guess of the parameters
param  = MPFIT('MYFUNCT', params,\\$
FUNCTARGS={X:x,Y:y,SIGMA_X:x_er,SIGMA_Y:y_er}, PERROR=p_er)```

The FUNCTARGS is the place where your data go into the function. And then you can find the results of the parameters from the fitting by print “param” and also the 1σ error from printing p_er. If you only have the error on “Y”, you can use a function like this:

```FUNCTION myfunct, p, x=x, y=y, er=er
a = p[0]
b = p[1]
c = p[2]
f = a*alog10(b)+c
RETURN, (y-f)/er
END```

And then using the same way described above to call “mpfit” and get the results and errors.

Here is also a PYTHON version of this procedure.