Received: from malur.postgresql.org ([217.196.149.56]) by arkaria.postgresql.org with esmtp (Exim 4.80) (envelope-from ) id 1XCBY9-0007QS-7y for pgsql-performance@arkaria.postgresql.org; Tue, 29 Jul 2014 17:50:57 +0000 Received: from localhost ([127.0.0.1] helo=postgresql.org) by malur.postgresql.org with smtp (Exim 4.80) (envelope-from ) id 1XCBY8-0006UX-Ex for pgsql-performance@arkaria.postgresql.org; Tue, 29 Jul 2014 17:50:56 +0000 Received: from magus.postgresql.org ([2a02:c0:301:0:ffff::29]) by malur.postgresql.org with esmtps (TLS1.2:DHE_RSA_AES_256_CBC_SHA256:256) (Exim 4.80) (envelope-from ) id 1XCBY6-0006UQ-Ow for pgsql-performance@postgresql.org; Tue, 29 Jul 2014 17:50:54 +0000 Received: from mail-wg0-x22e.google.com ([2a00:1450:400c:c00::22e]) by magus.postgresql.org with esmtps (TLS1.0:RSA_AES_256_CBC_SHA1:256) (Exim 4.80) (envelope-from ) id 1XCBY1-0001dZ-Mc for pgsql-performance@postgresql.org; Tue, 29 Jul 2014 17:50:53 +0000 Received: by mail-wg0-f46.google.com with SMTP id m15so5854wgh.17 for ; Tue, 29 Jul 2014 10:50:48 -0700 (PDT) DKIM-Signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=20120113; h=from:to:cc:references:in-reply-to:subject:date:message-id :mime-version:content-type:thread-index:content-language; bh=7GXsAUVmX93fh4yqgp1CHgg3UoGNRJ9VcUNVLQbiT4E=; b=lwvJw0lYLDDAZZGxRvPkpZfpHJ23I2t2TaITo80QvUkjO0zO5fHTWMYCcDUbDmB1fP /lecYar7RlvPYtO6jcMzhg4Bm904rAHb7s3Y08ntnZULC+9CSHSNudw0pvmOyIbBDzZg sQq0xmTmaAViArCXlF/g/zFN2c+GxusC3K8aZzXZXvLjrzAg4if76mHTQIZzFxctXhLv atKhRUh9OwGeYql/2GaBiQq0RHphzD8GAi+9x+wahSZqSbcRWczkqGxHpQ/K8GAf71Qi Zq+u1tHuW4kuN0LLeZ6caGN6rMJdUOH5G+ryZfxK3y3vrba/oEdNF7JvkGSuSUjbYeRT YUBg== X-Received: by 10.194.58.180 with SMTP id s20mr5278520wjq.119.1406656248246; Tue, 29 Jul 2014 10:50:48 -0700 (PDT) Received: from ssvnadvornik (ip-213-220-249-54.net.upcbroadband.cz. [213.220.249.54]) by mx.google.com with ESMTPSA id v4sm45732824wiz.16.2014.07.29.10.50.45 for (version=TLSv1 cipher=ECDHE-RSA-AES128-SHA bits=128/128); Tue, 29 Jul 2014 10:50:46 -0700 (PDT) From: =?UTF-8?B?SmnFmcOtIE7DoWR2b3Juw61r?= To: "'Sergey Karpov'" , "'Oleg Bartunov'" Cc: References: <01d201cfa8be$dd746880$985d3980$@gmail.com>

In-Reply-To: Subject: Re: Cursor + upsert (astronomical data) Date: Tue, 29 Jul 2014 19:50:43 +0200 Message-ID: <005001cfab55$a0ace910$e206bb30$@gmail.com> MIME-Version: 1.0 Content-Type: multipart/mixed; boundary="----=_NextPart_000_0051_01CFAB66.64454970" X-Mailer: Microsoft Outlook 15.0 Thread-Index: AQHSB3FUPUcQjIVcmXO3GnoQtcDS6gF+S2a+AU1F20ebm83EMA== Content-Language: cs X-Pg-Spam-Score: -2.0 (--) List-Archive: List-Help: List-ID: List-Owner: List-Post: List-Subscribe: List-Unsubscribe: X-Mailing-List: pgsql-performance Precedence: bulk Sender: pgsql-performance-owner@postgresql.org This is a multipart message in MIME format. ------=_NextPart_000_0051_01CFAB66.64454970 Content-Type: text/plain; charset="UTF-8" Content-Transfer-Encoding: quoted-printable Hello Sergey, Oh dear, I should have written you before several months :). I actually did= exactly what you suggest - it was actually Markus Demleitners (GAVO org.) = idea. The groupby in ipixcenter runs indeed less than an hour on the whole = table and is linear by definition so no problem with memory. I used I think= bit 11, which means when I proclaim that the ipix is a square and it has t= he number of steradians (one of your func. I don't remember which), it has = cca 0.7 arcsec side and cca 1.1 arcsec "half-diameter" which means it appro= ximately covers the circle of 1 arcsec radius with the center in the square= s center - if you understand what I mean. This approach and it's quite good results I actually presented on IVOA Inte= rop in Madrid this year (don't be afraid, I pointed out, that I use your Q3= C plugin for PostgreSQL).=20 The problem actually comes when you have more than 2 of these ipix "squares= " in line - when joining them, you don't know which ones to use. It causes = about 22% of these joined ipix centers to have collisions - when I take the= resulting ipix (joined from the ipix centers) and look back for all the or= iginal observations from which these ipix centers where constructed - they = don't belong to the same object. It is because when you count with the ipix= centers instead of the original observations - in most pathological case y= ou shift the observations from the bottom left corner to the top upper corn= er when joining with a square to your upper left. That means up to 2.2 arcs= ec shift which you can't backtrack. But. It definitely brought us closer to our goal - you don't have 3e8 obser= vations, you have cca 7e7 observations. It's still too much to process with= some quadratical algorithm, but it's something..=20 Maybe I could use it along with my *na=C3=AFve* sequential algorithm.. Coul= d you look at the script attached if you could see any obvious inefficient = statements? Thank you very much! You know even talking with someone about the problem b= rings you closer to the solution. P.S.: I am sorry if I got a little confused with your name when addresing y= ou - Are you the co-author of Q3C paper where the name is S. Koposov (+ O. = Bartunov)? Best Regards, Jiri Nadvornik -----Original Message----- From: Sergey Karpov [mailto:karpov.sv@gmail.com]=20 Sent: Monday, July 28, 2014 5:03 PM To: Oleg Bartunov Cc: Ji=C5=99=C3=AD N=C3=A1dvorn=C3=ADk; pgsql-performance@postgresql.org Subject: Re: [PERFORM] Cursor + upsert (astronomical data) Hi Jiri, I understand your problem (and I actually have exactly the same in my sky m= onitoring experiment). Unfortunately, I have no complete solution for it as= of now. I just may suggest you to look at q3c_ipixcenter() function from q3c which = returns the Q3C bin number for a given coordinates at a given binning lengt= h. Probably, you may first group your observations using this function, sel= ecting the binning depth to roughly match your desired angular resolution, = and then probably merge the (much shorter) list of resulting bin numbers (s= ort of common identifier you are looking for) for corner cases (when the sa= me cluster of observational points lies between two bins). I did not try th= is approach myself yet, but probably it will be significantly faster. Good luck! Sergey 2014-07-27 20:56 GMT+04:00 Oleg Bartunov : > Jiri, > > as I understand your problem is called crossmatch ? I attach pdf of=20 > our work in progress, where we compared several spatial indexing=20 > techniques, including postgis, q3c and new pgsphere. Sergey Karpov=20 > made these benchmarks. > > New pgsphere you may find here > https://github.com/akorotkov/pgsphere > > Oleg > > On Sat, Jul 26, 2014 at 1:46 PM, Ji=C5=99=C3=AD N=C3=A1dvorn=C3=ADk wrote: >> Hello guys. >> >> >> >> My issue kind of hits multiple topics, but the main question is about=20 >> performance. I think you need to understand the background a little=20 >> bit to be able to help me. So I will firstly define the problem and=20 >> my solutions to it and place the questions for you to the end of this me= ssage. >> >> >> >> Problem: >> >> I have a table of observations of objects on the sky. The most=20 >> important columns are the coordinates (x,y). All other columns in=20 >> there are just additional information about the observation. The=20 >> problem is that when I take an image of the same object on the sky=20 >> twice, the coordinates x,y won=E2=80=99t be the same, they will be only = close=20 >> to each other. My task is to generate a common identifier to all of=20 >> the observations of the same object and assign the observations to it=20 >> (N:1 relation). The criterium is that all of the observations which=20 >> are within 1 arcsec of this identifier are considered as the same=20 >> object. I keep the identifiers in a separate table (objcat) and have a f= oreign key in the observations table. >> >> The reason why I solve the performance issues here is that the table=20 >> of observations has atm cca 3e8 rows after 1.5 year of gathering the=20 >> data. The number growth is linear. >> >> >> >> Technical background: >> >> I=E2=80=99m trying to keep the whole algoritm in DB if possible because = I=20 >> have a good PostgreSQL plugin Q3C for indexing the coordinates of the=20 >> objects on the sky=20 >> (https://code.google.com/p/q3c/source/browse/README.q3c). It also has=20 >> quite a few stored procedures to look in that table for near=20 >> neighbors which is what I=E2=80=99m doing. The function q3c_join(x1, x2,= y1,=20 >> y2, radius) returns true if the object y is within radius of the=20 >> object x. It simply generates a list of index bitmap or queries with=20 >> the operators <=3D, >=3D which define the position on sky. Asking for th= e nearest neighbors are then only index scans. >> >> >> >> Solution: >> >> After lot of experimentation with clustering the objects and trying=20 >> to process them all together in one =E2=80=9Csilver-bullet=E2=80=9D SQL = query I=20 >> decided to use some more simple approach. The main problem with the=20 >> =E2=80=9Cprocess all at once approach=E2=80=9D is that the finding neigh= bours for=20 >> each observation is in definition quadratic and for 3e8 rows just=20 >> runs out of disk space (~TBs of memory for the temporary results). >> >> The simplest approach I could think of is that I process each row of=20 >> the 3e8 rows sequentially and ask: >> >> Do I have an identifier in the radius of 1 arcsec? >> >> No: Generate one and assign me to it. >> >> Yes: Update it and assigne me to it. The update is done as weighted=20 >> average =E2=80=93 I keep the number of how many observations the identif= ier=20 >> has been computed. The result is that the identifier will have=20 >> average coordinates of all the observations it identifies =E2=80=93 it w= ill be the center. >> >> >> >> So here I come with my procedure. It has 3 params. The first two are=20 >> range of oids to list in the table. Used for scaling and=20 >> parallelization of the algorithm. The third is the radius in which to se= arch for the neighbours. >> >> >> >> DROP TYPE IF EXISTS coords; >> >> CREATE TYPE coords AS ( >> >> raj2000 double precision, >> >> dej2000 double precision >> >> ); >> >> >> >> >> >> DROP FUNCTION IF EXISTS build_catalog(int,int,double precision); >> >> CREATE OR REPLACE FUNCTION build_catalog (fromOID int, toOID int,=20 >> radius double precision) >> >> RETURNS VOID AS $$ >> >> DECLARE >> >> cur1 CURSOR FOR >> >> SELECT >> >> raj2000, dej2000 >> >> FROM >> >> \schema.observation AS=20 >> obs >> >> WHERE >> >> obs.oid >=3D fromOID >> >> AND >> >> obs.oid < toOID; >> >> curr_raj2000 double precision; >> >> curr_dej2000 double precision; >> >> curr_coords_cat coords; >> >> cnt int; >> >> >> >> BEGIN >> >> /*SELECT current_setting('transaction_isolation') into tmp; >> >> raise notice 'Isolation level %', tmp;*/ >> >> OPEN cur1; >> >> cnt:=3D0; >> >> LOCK TABLE \schema.objcat IN SHARE ROW EXCLUSIVE MODE; >> >> LOOP >> >> FETCH cur1 INTO curr_raj2000, curr_dej2000; >> >> EXIT WHEN NOT found; >> >> >> >> WITH >> >> upsert >> >> AS >> >> (UPDATE >> >> \schema.objcat >> >> SET >> >> ipix_cat=3Dq3c_ang2ipix( >> >>=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20 >> (raj2000 * weight + curr_raj2000) / (weight + 1), >> >>=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20 >> (dej2000 * weight + curr_dej2000) / (weight + 1) >> >> ), >> >> raj2000 =3D (raj2000 *=20 >> weight + >> curr_raj2000) / (weight + 1), >> >> dej2000 =3D (dej2000 *=20 >> weight + >> curr_dej2000) / (weight + 1), >> >> weight=3Dweight+1 >> >> WHERE >> >> q3c_join(curr_raj2000,=20 >> curr_dej2000, >> >>=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20 >> raj2000, dej2000, >> >>=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20 >> radius) >> >> RETURNING *), >> >> ins AS >> >> (INSERT INTO >> >> \schema.objcat >> >> (ipix_cat, raj2000,=20 >> dej2000, >> weight) >> >> SELECT >> >>=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20= =20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20=20 >> (q3c_ang2ipix(curr_raj2000, curr_dej2000)), >> >> curr_raj2000, >> >> curr_dej2000, >> >> 1 >> >> WHERE NOT EXISTS >> >> (SELECT * FROM upsert) >> >> RETURNING *) >> >> UPDATE >> >> \schema.observation >> >> SET >> >> id_cat =3D (SELECT DISTINCT >> >> id_cat >> >> FROM >> >> upsert >> >> UNION >> >> SELECT >> >> id_cat >> >> FROM >> >> ins >> >> WHERE id_cat IS NOT=20 >> NULL >> >> LIMIT 1) >> >> WHERE CURRENT OF cur1; >> >> cnt:=3Dcnt+1; >> >> >> >> IF ((cnt % 100000 ) =3D 0) THEN >> >> RAISE NOTICE 'Processed % entries',=20 >> cnt; >> >> END IF; >> >> >> >> END LOOP; >> >> CLOSE cur1; >> >> END; >> >> $$ LANGUAGE plpgsql; >> >> >> >> Results: When I run the query only once (1 client thread), it runs=20 >> cca 1 mil rows per hour. Which is days for the whole dataset. When I=20 >> run it in parallel with that lock to ensure pessimistic=20 >> synchronization, it runs sequentially too J the other threads just=20 >> waiting. When I delete that lock and hope to solve the resulting=20 >> conflicts later, the ssd disk serves up to 4 threads relatively=20 >> effectively =E2=80=93 which can divide my days of time by 4 =E2=80=93 st= ill inacceptable. >> >> >> >> The reason is quite clear here =E2=80=93 I=E2=80=99m trying to write som= ething in one=20 >> cycle of the script to a table =E2=80=93 then in the following cycle I n= eed=20 >> to read that information. >> >> >> >> Questions for you: >> >> 1. The first question is if you can think of a better way how to do >> this and maybe if SQL is even capable of doing such thing =E2=80=93 or d= o I=20 >> have to do it in C? Would rewriting the SQL function to C help? >> >> 2. Could I somehow bend the commiting during the algorithm for my >> thing? Ensure that inside one cycle, the whole part of the=20 >> identifiers table would be kept in memory for faster lookups? >> >> 3. Why is this so slow? J It is comparable to the quadratic algori= thm >> in the terms of speed =E2=80=93 only does not use any memory. >> >> >> >> I tried to sum this up the best I could =E2=80=93 for more information p= lease=20 >> don=E2=80=99t hesitate to contact me. >> >> >> >> Thank you very much for even reading this far. >> >> >> >> Best Regards, >> >> >> >> Jiri Nadvornik >> >> >> >> ------=_NextPart_000_0051_01CFAB66.64454970 Content-Type: application/octet-stream; name="build_catalog_linear.sql" Content-Transfer-Encoding: quoted-printable Content-Disposition: attachment; filename="build_catalog_linear.sql" =EF=BB=BFalter table \schema.observation set with oids;=0A= alter table \schema.objcat drop column if exists id_cat cascade;=0A= alter table \schema.objcat add column id_cat serial;=0A= =0A= DROP TYPE IF EXISTS coords;=0A= CREATE TYPE coords AS (=0A= raj2000 double precision,=0A= dej2000 double precision=0A= );=0A= =0A= =0A= DROP FUNCTION IF EXISTS build_catalog(int,int,double precision);=0A= CREATE OR REPLACE FUNCTION build_catalog (fromOID int, toOID int, radius do= uble precision)=0A= RETURNS VOID AS $$=0A= DECLARE=20=0A= cur1 CURSOR FOR=20=0A= SELECT=20=0A= raj2000, dej2000=0A= FROM=0A= \schema.observation AS obs=0A= WHERE=20=0A= obs.oid >=3D fromOID=0A= AND=0A= obs.oid < toOID;=0A= curr_raj2000 double precision;=0A= curr_dej2000 double precision;=0A= curr_coords_cat coords;=0A= cnt int;=20=20=20=20=20=20=20=20=0A= =0A= BEGIN=20=0A= /*SELECT current_setting('transaction_isolation') into tmp;=0A= raise notice 'Isolation level %', tmp;*/=0A= OPEN cur1;=0A= cnt:=3D0;=0A= LOCK TABLE \schema.objcat IN SHARE ROW EXCLUSIVE MODE;=0A= LOOP=20=0A= FETCH cur1 INTO curr_raj2000, curr_dej2000;=0A= EXIT WHEN NOT found;=0A= =09=0A= WITH=20=0A= upsert=0A= AS=0A= (UPDATE=20=0A= \schema.objcat=0A= SET=20=0A= ipix_cat=3Dq3c_ang2ipix(=0A= (raj2000 * weight + curr_raj2000) / (weight + 1),=0A= (dej2000 * weight + curr_dej2000) / (weight + 1)=0A= ),=0A= raj2000 =3D (raj2000 * weight + curr_raj2000) / (weight + 1),=0A= dej2000 =3D (dej2000 * weight + curr_dej2000) / (weight + 1),=0A= weight=3Dweight+1=20=0A= WHERE=20=0A= q3c_join(curr_raj2000, curr_dej2000,=0A= raj2000, dej2000,=0A= radius)=0A= RETURNING *),=0A= ins AS=0A= (INSERT INTO=20=0A= \schema.objcat=20=0A= (ipix_cat, raj2000, dej2000, weight)=0A= SELECT=0A= (q3c_ang2ipix(curr_raj2000, curr_dej2000)),=0A= curr_raj2000,=0A= curr_dej2000,=0A= 1=0A= WHERE NOT EXISTS=0A= (SELECT * FROM upsert)=0A= RETURNING *)=0A= UPDATE=20=0A= \schema.observation=0A= SET=20=0A= id_cat =3D (SELECT DISTINCT=0A= id_cat=0A= FROM=0A= upsert=0A= UNION=20=0A= SELECT=20=0A= id_cat=0A= FROM=0A= ins=0A= WHERE id_cat IS NOT NULL=0A= LIMIT 1)=0A= WHERE CURRENT OF cur1;=0A= cnt:=3Dcnt+1;=0A= =0A= IF ((cnt % 100000 ) =3D 0) THEN=0A= RAISE NOTICE 'Processed % entries', cnt;=0A= END IF;=0A= =09=09=0A= END LOOP;=0A= CLOSE cur1;=0A= END;=0A= $$ LANGUAGE plpgsql;=0A= =0A= SELECT build_catalog(854536,880530, 2.77778e-4);= ------=_NextPart_000_0051_01CFAB66.64454970 Content-Type: text/plain Content-Disposition: inline Content-Transfer-Encoding: 8bit MIME-Version: 1.0 -- Sent via pgsql-performance mailing list (pgsql-performance@postgresql.org) To make changes to your subscription: http://www.postgresql.org/mailpref/pgsql-performance ------=_NextPart_000_0051_01CFAB66.64454970--