LAFs2ell

[yunyundong]

There were no description about functions LAFs2ell, LAFs2ellT, ells2LAFsT. Can you provide the reference about the relationship between LAF and ell? It is difficult to understand the relation only from the source code. Thank you in advance! @ @ducha-aiki

[ducha-aiki]

LAF is local affine frame, 2x2 matrix A (2x3, if including x and y)
It is multiplication of the rotation matrix (cos a sin a; -sin a cos a), scale (radius of the patch), anisotropic scale matrix (sc1 0; 0 1) and orientation of the anisotropic scaling axis (cos b sin b; -sin b cos b).
It is very convenient way of representing oriented elliptic features.

You can read about it here:
https://cw.fel.cvut.cz/b172/courses/mpv/labs/2_correspondence_problem/start

Ellipse is orientation-less representation of the ellipse by 3 coefficients: a, b, c.
It is standard output here http://www.robots.ox.ac.uk/~vgg/research/affine/detectors.html

Relationship between them is ellipse (a,b,c) = invsqrt(LAF)

T stand for pytorch function, no T means numpy function doing the same

[yunyundong]

Hi, the A23 returned by the function Ell2LAF(ell) is the multiplication of (sc1 0; 0 1) and (cos b sin b; -sin b cos b)?
I dont see the role of rectifyAffineTransformationUpIsUp() , can you give me a details?
Besides, if I have a local affine frame ,such as a=[1,2;3,4], using LAFs2ell function to convert it to the corresponding orientation-less ellipse ell_a, then I plugged the resulting ell_a to the function Ell2LAF , I can not get the same local affine frame a.

[a,b;b,c]=inv(sqrtm(LAF)) is right ? or LAF=inv(sqrtm([a,b;b,c]))

def Ell2LAF(ell):
A23 = np.zeros((2,3))
A23[0,2] = ell[0]
A23[1,2] = ell[1]
a = ell[2]
b = ell[3]
c = ell[4]
sc = np.sqrt(np.sqrt(ac - bb))
ia,ib,ic = invSqrt(a,b,c) #because sqrtm returns ::-1, ::-1 matrix, don`t know why
A = np.array([[ia, ib], [ib, ic]]) / sc
sc = np.sqrt(A[0,0] * A[1,1] - A[1,0] * A[0,1])
A23[0:2,0:2] = rectifyAffineTransformationUpIsUp(A / sc) * sc
return A23

In this function, I dont see the role of rectifyAffineTransformationUpIsUp. In fact, import scipy.linlag as la
, A =la.inv(la.sqrtm(ell))).

In this function LAFs2ell, its input LAF, is the

multiplication of the rotation matrix (cos a sin a; -sin a cos a), scale (radius of the patch), anisotropic scale matrix (sc1 0; 0 1) and orientation of the anisotropic scaling axis (cos b sin b; -sin b cos b)?

Thank you in advance. @ducha-aiki

[ducha-aiki]

A23 is not the result of multiplication, but it can be interpreted as such.
Operation laf-ell-laf is not invertible, because laf is 4 parameters, while ellipse is 3. So ellipse is unoriented ellipse, while LAF is oriented.
rectifyAffineTransformationUpIsUp cancels the orientation in LAF, making the orientation always equal to gravity vector, that is why it is called "up is up". After such rectification, LAF now effectively has also 3 parameters.

[yunyundong]

For a unit circle, it satisfy the condition x'x=1 (here, x is a column vector). Given an affine transformation matrix A=[a,b;c,d], applying affine transformation to x, so we have y=Ax. Now, y constitutes a ellipse. How can I solve the long axis and short axis of the ellipse? How can I solve the angle of ellipse? Are there materials to introduce the derivation of the solver? @ducha-aiki

[ducha-aiki]

Sorry, I don`t understand your question about solving and solvers, nor how is it related to affnet.
Take A = [a, b, c, d]. sqrt(det A) = sqrt(ac - bd) is isotropic scale.
Now A1 = A / det A. Take SVD (A1), and you will obtain (cos1, sin1; -sin1 cos1) * (t1 0; 0 t2) (cos2, sin2; -sin2 cos2).
t1, t2 = relative length of the axes of the ellipse, s.t. t1t2 = 1.
Is it what do you want?

[ducha-aiki]

@yunyundong See Section.2 here. Hope this helps
http://cmp.felk.cvut.cz/~chum/papers/perdoch-cvpr09.pdf

[yunyundong]

Take SVD (A1), and you will obtain (cos1, sin1; -sin1 cos1) * (t1 0; 0 t2) (cos2, sin2; -sin2 cos2).
t1, t2 = relative length of the axes of the ellipse, s.t. t1t2 = 1.

u,w,v=SVD(A1), according to the acquired u,w,v, how to calculate cos1,cos2,t1, t2?

[ducha-aiki]

u = (cos1, sin1; -sin1 cos1)
w = (t1 0; 0 t2)
v = (cos2, sin2; -sin2 cos2).

[yunyundong]

It is multiplication of the rotation matrix (cos a sin a; -sin a cos a), scale (radius of the patch), anisotropic scale matrix (sc1 0; 0 1) and orientation of the anisotropic scaling axis (cos b sin b; -sin b cos b).

so, here, cos1==cosa, cos2==cosb , are they right?

[ducha-aiki]

Yes. Or 1 = b, 2 = a. I don`t explicitly remember which rotation is which, but it is easy to figure out by trying them out.

[yunyundong]

It is multiplication of the rotation matrix (cos a sin a; -sin a cos a), scale (radius of the patch), anisotropic scale matrix (sc1 0; 0 1) and orientation of the anisotropic scaling axis (cos b sin b; -sin b cos b).

sc1 !=1, however t1*t2==1, what seems to be missing a scalar?

[yunyundong]

Take SVD (A1), and you will obtain (cos1, sin1; -sin1 cos1) * (t1 0; 0 t2) (cos2, sin2; -sin2 cos2).
t1, t2 = relative length of the axes of the ellipse, s.t. t1t2 = 1.
u = (cos1, sin1; -sin1 cos1)
w = (t1 0; 0 t2)
v = (cos2, sin2; -sin2 cos2).
where, det(A1)=1.

This is a very beautiful relationship. Where can I find the corresponding references. Thank you very much.

[ducha-aiki]

@yunyundong Missing scalar is in isotropic scaling scale, you can normalize anisotropic tilt matrix both ways. If you do it via SVD, as I described it in the last comments, it would be t1, t2, t1*t2 = 1.

However, it is not convenient to human to have (3 0; 0 1/3), because then it is not obvious from the first glance that one axis is 9 bigger than other. So better to have(9 0; 0 1), if you constructing such things from the elements, not by decomposition.

You can read about this from many sources. One of them is
http://www.ipol.im/pub/art/2011/my-asift/article_lr.pdf

[yunyundong]

so, sc1==t1/t2 (if t1>t2)

[ducha-aiki]

@yunyundong yes.
P.S. Could you please tell me for which task such details are important? Such things are derived through elementary matrix operations if needed, and lots of other parameterizations can be obtained having some task in mind. And giving the derivation nature of such things, I never actually remember exact order, signs, etc until the time of implementation, just general principle.

[yunyundong]

I am a newer to the representation of detected region. I only want to understand your source code thoroughly, but I failed. After your detailed interpretation, now I can understand your source code.
I want to use affnet to replace with sift algorithm to detect feature points in remote sensing images, is it suitable?
Besides, I want to know is there a detailed proof for a SVD of A ,where det(A)=1, the singular vector seemly is a rotation matrix.

[ducha-aiki]

@yunyundong It depends on what do you mean by "SIFT". Some times it understood as detector+orientation estimator + descriptor, as described in original paper. However, the detector from SIFT paper is difference of gaussians (DoG) and in computer vision community SIFT means only descriptor.

So, if you want to DETECT points, then you could use this example:
https://github.com/ducha-aiki/affnet/blob/master/examples/hesaffnet/hesaffnet.py

if instead you want to detect and DESCRIBE points and then match them with another image, I recommend you this example:
https://github.com/ducha-aiki/affnet/blob/master/examples/hesaffnet/WBS%20demo.ipynb

But such example is missing RANSAC part, so full matching pipeline with homography estimation is best from here:

https://github.com/ducha-aiki/mods-light-zmq

Regarding SVD, I don`t know where is the proof, but you can easily derive it by yourself, given that u and v are orthonormal matrices and all orthonormal matrices of 2x2 are rotation matrices.

P.S. I recommend you to read not only source code, but the paper itself, as well as Hessian-Affine paper. https://www.robots.ox.ac.uk/~vgg/research/affine/det_eval_files/mikolajczyk_ijcv2004.pdf

[yunyundong]

@ducha-aiki Thank you very much for your kind help. I will read them. Thank you again.

[yunyundong]

def LAFs2ell(in_LAFs):
LAFs = convertLAFs_to_A23format(in_LAFs)
ellipses = np.zeros((len(LAFs),5))
for i in range(len(LAFs)):
LAF = deepcopy(LAFs[i,:,:])
scale = np.sqrt(LAF[0,0]*LAF[1,1] - LAF[0,1]*LAF[1, 0] + 1e-10)
u, W, v = np.linalg.svd(LAF[0:2,0:2] / scale, full_matrices=True)
W[0] = 1. / (W[0]*W[0]scalescale)
W[1] = 1. / (W[1]*W[1]scalescale)
A = np.matmul(np.matmul(u, np.diag(W)), u.transpose())
ellipses[i,0] = LAF[0,2]
ellipses[i,1] = LAF[1,2]
ellipses[i,2] = A[0,0]
ellipses[i,3] = A[0,1]
ellipses[i,4] = A[1,1]
return ellipses

@ducha-aiki I want to know is there some special tricks for the function LAFs2ell?
LAF denotes the Affine transform A,, the attained ell is S=[s11,s12;s12,s22]. They should satisfy the relationshp AA'=S . So, why you do not use S=AA' directly?

import scipy.linlag as la
A=np.array([[9,1],[1,9]]) # denotes an affine transform
eA=LAFs2ell(A) #
u,w,b=la.svd(eA)
eA2=np.matmul(A,A.transpose()) # eA2=AA'
u1,w1,v1=lat.svd(eA2)

I found w[0]/w[1]==w1[0]/w1[1], so this shows that eA and eA2 denote the same shape of ellipse, however, eA/eA2!=scalar*[[1,1],[1,1]]. where i am wrong?

[ducha-aiki]

Good point, may be I could use such transfer directly. One important thing that I do not need ellipse format anywhere except for benchmark purposes, as benchmark scripts expect such.
It is the only point, where I need a,b,c format, so I have no incentive to optimize such thing at all. For the rest purposes I have, I work with LAFs directly.

[ducha-aiki]

I think, that your mistake is that SVD(matrix) != scalar * (SVD(matrix squared))

[yunyundong]

Hi, @ducha-aiki
What's the shape of nmsed_resp in the function https://github.com/ducha-aiki/affnet/blob/f46f0dcee547fb571c1c2d20adcd5e85cd6317b0/HandCraftedModules.py#L260 ? nmsed_resp, idxs = torch.topk(nmsed_resp, k = num_features, dim = 0); the shape of input nmsed_resp should be bachsize*3*height*width. So, the shape of output nmsed_resp is batchsize*3*height*k ? Is it right?
What's the role of the small piece of code? https://github.com/ducha-aiki/affnet/blob/f46f0dcee547fb571c1c2d20adcd5e85cd6317b0/HandCraftedModules.py#L278

Especially for the line sc_y_x = F.conv2d(resp3d, self.grid, padding = 1) / (F.conv2d(resp3d, self.grid_ones, padding = 1) + 1e-8)

[ducha-aiki]

@yunyundong answered in #15
If your new question is related not to the current issue - open new issue, please.